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FOREWORD 


This  report  was  prepared  by  the  Environmental  Control  Branch  of  the  Booster  and 
Power  Division.  Dyna-Soar  Engineering  Office,  Deputy  for  Engineering,  Aeronautical 
Systems  Division.  Wright-Patterson  Air  Force  Base,  Ohio,  under  Weapon  System  620A, 

The  work  was  performed  for  the  University  of  Buffalo,  Buffalo,  New  York,  in  partial 
fulfillment  of  the  requirement  for  the  degree  of  Master  of  Science  in  the  field  of 
Mechanical  Engineering.  Professor  Howard  E.  Strauss  served  as  graduate  advisor.  The 
University  of  Buffalo  has  granted  permission  to  the  USAF  to  publish  and  distribute  the 
thesis  as  an  ASD  Technical  Documentary  Report.  Only  those  changes  necessary  to  make 
the  thesis  meet  the  requirements  of  an  ASD  Technical  Documentary  Report  have  been 
made. 
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ABSTRACT 


Atmospheric  re-entry  of  earth -orbital,  hypersonic  glide  vehicles  creates  thermal 
problems.  The  heat  affects  not  only  the  materials  and  construction  of  the  airframe  but 
also  the  crew  and  various  subsystems  of  the  vehicle.  Successful  solution  of  these  problems 
depends  upon  the  development  of  an  effective  thermal  protective  concept,  which  will  also 
give  the  designer  some  latitude  in  his  design  philosophy.  The  role  of  the  protective 
system  is  to  slgnificrntly  attenuate  the  influx  of  heat  that  is  aerodynamically  generated 
within  the  surrounding  boundary  layer.  Attenuation  is  accomplished  by  combining  external 
radiation  shielding  elements  with  backup  insulation  materials  and  an  appropriate  cooling 
system. 

Analytical  procedures  are  presented  for  determining  significant  system  parameters 
by  transforming  the  differential  heat  conduction  or  diffusion  equation  into  an  algebraic 
expression  by  employing  the  calculus  of  finite  differences.  The  adaptation  of  the  resulting 
equation  to  digital  computer  programming  is  discussed,  and  numerical  results  are 
presented  to  indicate  systems  of  minimum  weight. 
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1.0  INTRODUCTION 

1.1  General 

Although  aeronautical  technology  offers  a  number  of  techniques  for  the  successful 
recovery  of  rocket -boosted  orbital  vehicles,  the  most  attractive  philosophy  is  derived 
from  the  principles  of  lifting  body  re-entry.  Unlike  trajectories  characteristic  of 
ballistic  or  semiballistic  shapes  that  are  controllable  only  to  the  extent  of  orbit  ejection 
sequencing,  a  lifting  body  or  glide  vehicle  can  be  piloted  to  a  preselected  site  where  a 
conventional  landing  is  made.  Atmospheric  re-entry  descent,  within  limited  parameters, 
depends  only  on  the  judicious  use  of  the  kinetic  and  potential  energy  of  the  vehicle. 

The  maneuvering  capability  of  this  type  of  vehicle  during  descent  can  be  defined  within 
velocity-altitude  limitations  of  a  flight  corridor  that  is  bounded  by  aerodynamic,  aero- 
elastic,  structural,  and  heating  considerations.  The  upper  limit  is  a  result  of  instability 
of  the  vehicle  at  high  angles  of  attack  or  a  maximum  attainable  aerodynamic  lift,  while 
the  lower  limit  is  a  function  of  maximum  heating  rates,  allowable  acceleration  loads,  or 
permissible  dynamic  pressure.  An  example  of  a  glide-vehicle  flight  corridor  is  shown  in 
Figure  1. 


Figure  1.  Flight  Corridor 
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Various  thermal  protection  concepts  affect  both  of  these  limits  for  the  following 
reasons: 

a.  Aerodynamic  lift  is  affected  by  wing  loading  or  by  the  weight  of  the  thermal 
protection  system. 

b.  Load  carrying  capability  is  a  function  of  structural  temperature,  which  is 
influenced  by  the  effectiveness  of  the  thermal  protection  system. 

As  a  result,  thermal  protection  system  optimization,  that  is,  adequate  protection  for 
minimum  weight,  contributes  to  a  broader  altitude  and  velocity  operational  capability. 

Glide  bodies,  unlike  ballistic  shapes,  experience  long  re-entry  times  while  exposed  to 
relatively  low  aerodynamic-heating  rates.  This  he,at  is  transferred  by  boundary-layer 
convection  to  the  surface  of  the  vehicle  where  it  can  be  radiated  to  space,  absorbed  by 
the  outer  skin,  and  transmitted  to  internal  sections.  Considerable  advantage  might  be 
realized  in  external  dissipation  of  this  convective  heat,  since  the  technique  does  not 
result  in  additional  vehicle  weight,  but  is  a  function  merely  of  the  emissive  power  of  the 
outside  surface.  A  knowledge  of  this  factor  along  with  the  rate  at  which  heat  is  aero- 
dynamically  generated  will  enable,  by  use  of  the  Stefan-Boltzmann  relationship,  an 
approximation  of  the  equilibrium  temperature  of  the  exterior  of  the  vehicle. 

Since  the  surface  equilibrium  temperatures  during  re-entry  flight  will  exceed  the 
internal  conditioned  compartment  temperature,  heat  will  penetrate  to  areas  of  compart¬ 
ment  location.  The  modes  of  transfer,  which  include  solid  conduction,  radiation,  gaseous 
conduction,  and  convection,  must  be  minimized  to  a  degree  predicated  by  the  relative 
contribution  of  each  to  the  total  heat  influx. 

Heat  conduction  that  passes  energy  directly  through  the  sections  of  a  solid  in  contact 
with  the  hotter  surfaces  can  be  lowered  by  using  materials  that  are  classified  as  poor 
thermal  conductors  and  by  reducing  the  cross-sectional  area  normal  to  the  flow  path. 

The  latter  may  be  achieved  by  incorporation  of  a  large  number  of  gaseous  spaces,  since 
conduction  through  a  gas  is  much  less  than  conduction  through  a  solid.  Another  technique 
would  involve  breaking  the  solid  down  into  a  collection  of  fine  particles,  which,  by  forming 
many  points  of  contact  within  the  material,  would  result  in  an  increased  resistance  to  the 
flow  of  heat. 

Heat  transfer  by  radiation  from  the  surface  of  a  hot  body  occurs  through  a  mechanism 
similar  to  electromagnetic  wave  phenomena  and  is,  in  fact,  frequency  dependent  and 
intimately  connected  to  that  of  light  transmission.  Geometric  orientation  of  the  surfaces 
relative  to  the  incident  radiation  and  techniques  for  reducing  surface  emissivity  both 
contribute  to  the  attenuation  of  this  process. 

The  transfer  of  heat  by  conduction  within  a  gas  is  described  as  diffusion  process  in 
which  the  molecules  are  in  motion  between  warm  and  cool  areas.  In  addition,  an 
exchange  in  kinetic  energy  occurs  as  wandering  molecules  collide  with  one  another. 

Gases  in  a  free  unconflned  state  exhibit  virtually  no  change  in  thermal  conductivity  as 
a  result  of  pressure  changes,  because  their  average  path  length  decreases  as  the 
number  of  molecules  increase  in  proportion  to  the  pressure.  On  the  other  hand,  confined 
gases  show  lower  thermal  conductivities  at  reduced  pressures  because  the  number  of 
molecules  available  for  transporting  heat  energy  approaches  zero.  This  principle  can  be 
applied  to  reduce  the  effects  of  gas  conduction  by  using  the  ambient  pressures  existing 
at  altitude  during  re-entry. 
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Heat  convection  takes  place  as  a  result  of  a  free -molecular  mixing  motion  initiated  by 
temperature  differences  within  a  gas.  If  the  gas  is  contained  within  small  spaces,  the 
motion  of  molecules  will  be  confined  to  an  extent  negating  heat  convection. 

If  during  re-entry  the  temperatures  within  the  pressurized  areas  of  the  vehicle  are 
maintained  essentially  at  equilibrium  without  using  compartment  conditioning  equipment 
to  absorb  heat  penetrating  from  the  external  surface,  then  this  function  must  be  assigned 
to  the  thermal  protection  system.  The  feasibility  of  this  direction  is  justified  in  view  of 
available  techniques  that  are  more  effective  than  those  that  depend  on  the  circulation  of 
a  gaseous  cooling  media  within  a  relatively  large  volume.  A  description  ofthe  aspects  of 
thermal  protection  systems  is  presented  in  the  following  section. 

1.2  Thermal  Protection  Concepts 

1.2.1  General 


Systems  protecting  the  internal  structure  or  compartments  from  the  effects  of  aero¬ 
dynamic  heating  consist  of  three  major  integrated  elements.  First,  the  exterior  surface 
forming  high  emissivity  heat-shield  segments  transmits  to  the  surrounding  environment 
a  great  deal  of  the  heat  reaching  the  vehicle  from  the  boundary  layer.  The  effectiveness 
of  this  scheme  in  minimizing  internal  heat  penetration  is  shown  in  Figure  2.  The  second 
element,  a  lightweight  thermal  insulation,  attenuates  the  remaining  heat  conducted  to 
areas  surrounding  the  environmentally  controlled  compartments.  At  this  location,  the 
heat  is  removed  by  a  cooling  system,  which  serves  as  the  remaining  element. 


.AVERAGE  SURFACE  TEMPERATURE  (DEGREES  FXIO"*) 


Figure  2.  Percent  of  Total  Heat  Removed  by  Cooling  System 
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1.2.2  Radiation  and  Ablation  Cooling 

The  primary  factors  associated  with  the  development  of  radiation  heat  shields  are 
those  governing  the  selection  of  materials,  emissivity  of  the  external  surface,  and  the 
evolution  of  a  suitable  design  configuration.  Criteria  forming  a  basis  for  material 
selection  include  oxidation  resistance  and  strength  requirements  that  are  compatible 
with  the  specified  service  life  and  the  external  environments.  As  a  result  of  the  latter, 
surface -temperature  histories  of  the  vehicle  are  predictable  functions  of  instantaneous 
heat  generation  rates  during  re-entry  and  the  emissivity  of  the  surface.  The  relation¬ 
ship  between  heating  rates,  emissivity,  and  heat-shield  equilibrium  temperature  is 
shown  in  Figure  3.  Details  of  design  configurations  are  concerned  more  with  the 
incorporation  of  stiffness  and  thermal  stress  alleviating  features  into  lightweight  heat 
shields  to  preclude  panel  buckling. 


Figure  3.  Effect  of  Surface  Emissivity  and  Heat  Transfer  Rate  on  Equilibrium  Temperature 


Vehicle  surfaces  that  can  withstand  temperatures  up  to  2000'’F  can  be  constructed 
with  super  alloy  materials  using  current  manufacturing  methods  and  materials  technology. 
As  metallurgical  research  uncovers  new  materials,  this  limit  might  be  extended  to 
3000°F.  Refractory  metals,  for  example,  are  generally  useful  to  2700°F  if  effective 
oxidation-resistant  coatings  are  available  and  high  surface  emissivities  can  be  main¬ 
tained.  Furthermore,  refractory  metals  are  not  suitable  for  many  structural  applica¬ 
tions;  hence,  aerodynamic  loads  must  be  transmitted  to  more  ductile  members.  Other 
refractory  materials  such  as  ceramics  and  graphites  possess  excellent  elevated  tempera¬ 
ture  characteristics  which,  when  integrated  with  suitable  design  concepts,  also  can  be 
applied  to  exterior  surface  panels. 
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Still  another  approach  is  that  of  using  an  ablative  surface  that  is  combined  with 
radiation  shields  to  give  a  re-entry  vehicle  capability  for  returning  from  superorbital 
missions.  Under  these  circumstances,  trajectories  that  result  in  heating  rates  that  are 
associated  with  the  ablation  process  must  be  considered. 

1.2.3  Thermal  Insulations  and  Cooling  Systems 


1.2. 3.1  General 

Radiation  cooling  that  provides  a  high  surface  emissivity  for  a  low  weight,  and 
materials  characterized  by  relatively  high  heats  of  ablation,  although  effective,  do  not 
completely  solve  the  problems  that  are  associated  with  thermal  protection  for  re-entry 
heating.  As  a  result  of  extended  flight  through  the  earth’s  atmosphere,  temperature 
gradients  are  established  that  cause  heat  to  penetrate  the  external  surface.  This  heat 
must  be  dissipated  to  maintain  environments  that  are  compatible  with  those  required 
for  human  comfort  and  equipment  operation. 

Several  approaches  can  be  considered  for  removing  this  heat  or  confining  it  to  a 
region  between  the  compartment  shell  and  the  outer  surface  of  the  vehicle. 

The  methods  that  involve  either  the  principles  of  heat  storage  or  the  vaporization  of  a 
liquid  are  as  follows; 

a.  Thermal  Insulation  systems 

b.  Liquid  heat  sink  systems 

c.  Systems  that  combine  thermal  insulation  and  liquid  heat  sinks. 

Low  thermai  diffusivity  insuiations  can  delay  the  influx  of  aerodynamic  hear  until 
completion  of  the  mission.  Two  obvious  physical-property  requirements  for  these 
materials  are  a  low  thermal  conductivity  and  high  specific  heat  at  the  temperature  and 
pressure  levels  that  result  from  re-entry  trajectories.  Insulations  available  from 
industry,  which  can  retain  their  form  at  elevated  temperatures  and  low  pressures  do 
not  display  values  for  the  properties  just  mentioned  to  permit  their  effective  application. 
Denser  forms  and  thicker  sections  of  the  materials  must  be  employed  as  compensating 
measures  to  obtain  the  desired  temperature-time  response  throughout  the  system. 

Protection  systems  using  the  latent  heat  of  vaporization  of  suitable  coolants  represent 
a  technological  approach  for  removing  boundary-layer  heat  that  penetrates  the  exterior 
surface.  The  total  heat  removed  or  the  weight  of  the  coolant  expanded  during  a  re-entry 
trajectory  is  a  function  of  the  thermal  conductance  between  the  heated  and  cooled  surfaces. 
A  convenient  means  for  reducing  the  conductance  involves  the  use  of  insulating  materials 
and  form  a  protection  system  configuration  when  combined  with  a  cooling  system. 

Conformation  of  the  schemes  selected  for  further  review  and  optimization  resulted 
from  the  use  of  relationships  derived  in  Reference  4  of  the  Bibliography.  Solutions  are 
presented  in  Figure  4  as  a  result  of  considering  step -function  temperature  inputs  that 
are  equivalent  to  mean  surface  temperatures  associated  with  re-entry  flight.  Weights 
for  combined  systems  were  determined  from  a  steady -state  heat  balance  that  equated 
the  heat  conducted  through  the  Insulation  to  that  removed  by  the  coolant.  The  results 
verify  that  the  most  effective  thermal -protection  systems  for  lifting  body  re-entry  vehicles 
combine  insulations  with  expendable  coolants. 
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1.2. 3. 2  Insulation  Materials 


The  selection  of  a  thermal  insulation  involves  a  comparison  of  properties  that  indicates 
the  potential  capabilities  of  the  materials.  The  insulating  materials  must  provide  not  only 
a  maximum  resistance  to  the  flow  of  heat  but  they  must  be  light  in  weight  and  demonstrate 
the  ability  to  endure  mission  environments  with  no  evidence  of  degradation  that  is  detri¬ 
mental  to  service  life. 

A  measure  of  the  relative  effectiveness  of  insulations  may  be  derived  by  expressing 
the  weight  of  the  material  as  a  function  that  depends  on  the  physical  properties.  If 
effectiveness  can  be  defined  as  the  reciprocal  of  the  weight  of  insulation,  the  function 
will  be  equal  to  the  inverse  of  the  density-thermal  conductivity  product.  The  effects  of 
altitude  on  this  product  should  be  considered  when  thermal  conductivity  values  are 
selected. 

A  comparison  of  several  commercially  available,  uncontained  insulating  materials  is 
shown  in  Figure  5.  The  possibility  exists  that  material  characteristics  dictate  contain¬ 
ment  concepts  that  significantly  affect  the  thermal  conductance  of  the  system.  For  these 
cases,  an  apparent  thermal  conductivity  that  includes  containment  effects  should  be  used 
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Figure  5.  Comparison  of  the  Effectiveness  of  hisulations 
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for  a  more  accurate  comparison.  These  effects  are  included  simply  by  a  weighted  adjust¬ 
ment  of  both  the  thermal  conductivity  and  density. 

Additional  important  characteristics  of  insulating  materials  are  the  following: 

a.  Adequate  compressive  strength  to  insure  retention  and  prevent  settling; 

b.  Softening  point  above  the  maximum  operating  temperature  to  prevent  a  loss  of 
compressive  strength  and  particle  shape; 

c.  Minimum  contact  area  between  individual  particles  to  reduce  solid  conduction; 

d.  Absence  of  sintering  at  operating  temperatures  to  prevent  degradation  of 
material  properties; 

e.  Ability  to  attenuate  radiant  heat  transfer  by  particle  orientation  and  emissive 
characteristics,  which  cause  sufficient  radiation  blockage; 

f.  Particle  size  distribution  and  construction  to  minimize  compaction  under  vibra¬ 
tional  or  gravitational  forces; 

g.  Low  vapor  pressure  at  maximum  operating  temperatures  to  preclude  loss  of  the 
material; 

h.  Particle  size  that  is  so  small  that  the  mean  free  path  of  the  air  molecules  be¬ 
comes  larger  than  the  particle  spacing  at  moderate  vacuums. 

Powder  materials  and  fibrous  forms  of  insulation  have  physical  characteristics  that 
effectively  reduce  each  of  the  contributing  modes  of  heat  transfer  and  satisfy  the 
preceding  requirements.  Gas  conduction,  for  example,  is  lowered  by  the  orientation  of 
fine  particles  to  create  pore  sizes,  which,  at  moderate  vacuums,  are  smaller  than  the 
mean  free  path  that  the  gas  molecules  travel  before  colliding.  Constituent  materials  for 
the  insulation  are  chosen  on  the  basis  of  their  ability  to  attenuate  heat  flow  through  the 
materials.  Major  constituents  such  as  alumina,  potassium  titanate,  and  fibrous  asbestos 
have  excellent  insulating  qualities,  which,  when  coupled  with  the  many  joint  resistances  at 
the  points  where  particles  contact  one  another,  resist  the  flow  of  heat  by  solid  conduction. 
Thermal  radiation  can  be  reduced  by  the  addition  of  both  radiation  absorbing  particles 
(such  as  zirconia,  carbon,  and  silicon  nitride)  and  radiation  reflecting  or  scattering 
particles  (such  as  aluminum,  platinum,  or  tantalum  flakes).  The  effect  of  temperature 
and  altitude  on  the  thermal  conductivity  of  an  insulation  powder  is  shown  in  Figure  6. 

Although  fine  powders  permit  a  wider  flexibility  in  selecting  ingredients  for  specific 
temperature  applications  than  do  fibrous  mats,  problems  that  are  associated  with  the 
packaging  of  these  materials  tend  to  overshadow  any  relative  advantages.  Specific 
problems  include  the  following: 

a.  Compaction  and  settling  of  constituents  under  vibratory  loads, 

b.  Excessive  resistance  to  air  flow  during  boost  venting, 

c.  Development  of  filter  configurations  to  retain  fine  particles, 

d.  Insulation  package  buckling  due  to  thermal  stress. 
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The  more  common  fibrous  insulations  contain  either  aluminum  silicate  or  silica  fibers, 
which  are  useful  to  2000°F  These  are  commercially  available  in  various  forms,  namely, 
bulk,  batt,  blanket,  cloth,  molded  shapes,  and  paper.  The  diameters  of  the  fibers  are  in  the 
range  from  0.00003  to  0.0004  inch  and  the  densities  range  from  3  to  26  pounds  per  square 
foot.  The  thermal  conductivity  for  one  of  the  fibrous  quartz  insulations  is  shown  in 
F  igure  7  as  a  function  of  temperature  and  altitude. 


Figure  6.  Effect  of  Temperature  and  Altitude  on  Thermal  Conductivity  of  Powdered  Mixture 
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Figure  7.  Effect  of  Temperature  and  Altitude  on  Thermal  Conductivity  of  Fibrous  Quartz  Insulation 
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1.2. 3. 3  Cooling  Systems 

Three  cooling  system  concepts  appear  feasible  for  use  in  removing  the  aerodynamic 
heat  that  is  transm.itted  through  the  insulation.  These  systems  are  defined  in  the 
succeeding  paragraphs. 

1.2. 3. 3.1  Open-Cycle  Systems 

The  coolant  of  the  open-cycle  systems  (a  99-percent  water,  1 -percent  gel  mixture)  is 
contained  by  a  sponge  material  within  passages  that  are  attached  to  a  thin  aluminum 
sheet.  The  assembled  panels  are  mounted  to  the  pilot  and  payload  compartments;  they 
enclose  areas  containing  equipment  that  cannot  tolerate  exposure  to  the  hot  surface 
panels.  When  transmitted  through  the  protective  insulation,  heat  is  absorbed  by  vaporizing 
the  coolant.  The  steam  generated  is  transported  to  an  overboard  dump  where  it  is 
exhausted  through  a  control  valve,  which  functions  primarily  to  regulate  the  coolant 
temperature.  Another  advantage  of  a  pressurized  system  is  the  elimination  of  coolant 
vaporization  during  orbit  since  the  vapor  pressure  of  the  gel  would  exceed  the  pressure 
within  the  passages.  A  schematic  of  this  system  is  shown  in  Figure  8. 


1.2. 3. 3. 2  Closed-Cycle  Systems 


The  compartment  shell  or  primary  structure  of  the  closed-cycle  systems,  constructed 
from  an  aluminum  alloy,  contains  an  integral  tube  circuit  that  forms  a  closed  transport 
loop  through  which  a  coolant  is  circulated.  During  re-entry,  heat  transmitted  through 
the  insulation  is  absorbed  by  the  circulating  fluid  and  transferred  in  a  remote  heat 
exchanger  to  an  expendable  heat  sink  that  uses  its  heat  of  vaporization  in  the  process. 
The  rate  of  flow  of  the  expendable  coolant  through  the  heat  exchanger  regulates  the 
temperature  in  the  closed  loop.  Figure  9  shows  a  schematic  of  a  closed-cycle  system. 
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Figure  9.  Closed-Cycle  Cooling  System 


1.2. 3. 3. 3  Combined  Open-Cycle  and  Closed-Cycle  Systems 

The  combined  open-cycle  and  closed-cycle  systems  can  be  used  advantageously  to 
enhance  the  overall  configuration  of  the  system.  In  general,  the  closed-cycle  systems 
can  be  adapted  readily  to  large  surface  areas  and  thermally  complex  sections  where 
local  conductance  values  vary  significantly.  Areas  where  appreciable  heat  conduction 
results  from  structural  framework  that  support  hatches  and  windshields,  and  members 
that  extend  between  the  outer  heat  shields  and  the  cooled  inner  structure  or  compartment 
walls  are  examples  of  this  case.  On  the  other  hand,  remote  equipment  and  movable 
surfaces  might  be  more  effectively  cooled  by  an  open  system  to  eliminate  long  or 
flexible  connections,  which  would  be  more  susceptible  to  fatigue. 
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1.3  System  Optimization 

In  the  final  analysis,  many  design  parameters  might  influence  optimization  procedures. 
These  include  allowable  material  temperatures,  reliability  considerations,  volume  limita¬ 
tions,  flexibility  and  fabrication  requirements  of  the  system,  systems  that  are  readily 
accessible  for  servicing,  andeaseof  adaptation  to  straightforward  preflight  checkout 
procedures.  In  addition,  the  thermal  protection  system  must  be  able  to  survive  environ¬ 
ments  that  are  associated  with  ground  support  as  well  as  all  phases  of  an  orbital  mission. 

The  integration  of  the  components  of  thermal  protection  systems  into  a  system  of 
minimum  weight  can  be  performed  after  the  type  of  insulation  and  cooling  system  have 
been  selected.  Essentially,  the  optimization  process  involves  a  tradeoff  between  insula¬ 
tion  and  coolant  requirements  to  arrive  at  a  combination  of  minimum  weight.  The 
following  sections  present  the  analytical  techniques  associated  with  this  process. 
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2.0  SUMMARY 


A  thermal  analysis  and  an  analytical  procedure  of  practical  significance  are  presented 
and  applied  to  establish  meaningful  design  parameters  for  thermally  protected  glide 
vehicles  that  re-enter  the  atmosphere.  With  only  slight  modification,  as  dictated  by 
particular  boundary  conditions,  the  techniques  developed  may  be  useful  for  obtaining 
solutions  to  a  number  of  problems  in  the  field  of  heat  transfer  where  temperatures  are 
functions  of  spatial  and  time  coordinates. 

The  general  differential  heat  conduction  or  diffusion  equation,  extended  to  include 
temperature -dependent  thermal  properties,  was  transformed  into  an  approximate 
algebraic  expression  by  employing  the  calculus  of  finite  differences.  The  difference 
terms  repiacing  the  partiai  derivatives  were  the  first  approximations  of  the  derivatives 
of  polynomials,  obtained  from  the  Gregory-Newton  and  Stirling  interpolation  formulas, 
which  represented  the  temperature  at  any  point  in  the  coordinate  system. 

Three  realistic  thermal  protection  systems  were  analyzed  by  the  application  of  this 
general  difference  equation  together  with  the  appropriate  boundary  conditions.  Three 
configurations  were  selected:  an  insulated  and  cooled  structure,  an  insulated  and  cooled 
compartment  combined  '.vith  a  “hot”  structure,  and  an  internally  radiation-cooled 
structure.  The  boundary  value  problems  in  each  case  were  reduced  to  expressions  that 
enabled  a  prediction  of  the  temperature  history  throughout  the  configuration,  instantaneous 
rates  of  heat  transfer  to  the  cooling  system,  as  well  as  the  total  coolant  requirements  that 
were  based  on  the  amount  of  heat  absorbed  during  a  typical  descent  trajectory. 

The  equations  resulting  from  the  analysis  of  the  insulated  and  cooled  structure  were 
programmed  on  an  IBM  7090  computer  and  numerical  solutions  were  obtained.  Tempera¬ 
ture  surveys  for  this  case  are  presented  in  the  results  and  the  significance  of  optimizing 
the  design  is  duscussed. 
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3.0  THERMAL  ANALYSIS 


3.1  Problem  Statement 

Optimization  procedures  establish,  by  an  orderly  analytical  process,  the  insulation  and 
cooling  requirements  for  thermal  protection  systems  of  minimum  weight.  These  require¬ 
ments  are  set  forth  by  arranging  the  problem  solutions  in  an  expression  for  the  total 
system  weight  as  a  function  of  insulation  thickness  and  then  by  selecting  the  thickness, 
which  corresponds  to  the  minimum  combined  weights.  The  techniques  involved  in 
generating  sufficient  data  to  adequately  define  this  function  will  be  presented  in  this 
section. 

The  total  weight  of  the  system  can  be  divided  into  two  basic  elements;  the  weight  of  the 
insulation  that  is  assumed  as  an  input  to  the  problem  and  the  coolant  weight  that  is 
expended  as  a  result  of  absorbing  heat,  which  is  conducted  through  the  insulation  from 
the  surface  of  the  vehicle.  The  latter  quantity  can  be  determined  from  the  solution  of  the 
boundary  value  problem  for  conditions  that  define  the  time-dependent  flow  of  heat 
through  the  insulation.  One  condition  on  U,  which  is  a  familiar  parabolic  partial  differential 
equation,  can  be  written  in  the  form 


U,  =  G  (U)  (1) 

where  the  right-hand  term  is  a  second  order  elliptic  partial  differential  operator  with 
either  one,  two,  or  three  independent  space  variables.  Now  if  R  is  the  region  of  the  space 
variable  or  variables  bounded  by  D,  then  the  values  of  U  in  R  can  be  determined  for  all 
t  >  to  after  U  in  R  is  specified  for  t  =  to  and  the  values  of  U  on  D  are  known  for  all  t  >to. 

Since  obtaining  solutions  to  this  problem  by  analytical  methods  is  not  practical,  an 
approximate  numerical  method  can  be  used  which  replaces  the  partial  differential  equation 
by  a  partial  difference  equation.  This  method,  described  here,  is  the  method  of  finite 
differences. 

A  grid  network  of  regularly  spaced  straight  lines,  each  parallel  to  one  of  the  coordinate 
axis,  is  superimposed  over  the  region  R.  In  addition,  a  time  increment,  p,  is  introduced, 
which  essentially  represents  a  new  dimension  to  the  mesh.  Now  instead  of  solving  a 
complex  differential  equation  for  all  values  in  a  continuous  region,  only  approximate 
values  of  the  solution  must  be  obtained  at  the  mesh  points  that  are  formed  by  the 
intersections  of  the  network  lines  with  each  other  and  with  the  boundary  of  the 
region. 

The  method  of  obtaining  the  difference  equation  for  each  of  the  interior  mesh  points 
involves  a  transformation  of  the  differential  equation  by  replacing  the  partial  derivatives 
with  corresponding  difference  quotients.  As  a  result,  the  problem  can  be  formulated 
into  a  series  of  linear  algebraic  equations  where  each  defines  the  conditions  at  a 
particular  interior  mesh  point.  In  the  case  of  the  boundary  value  problems  that  are 
considered  in  the  following  sections,  the  solutions  will  involve  a  system  of  n  linear 
algebraic  equations  with  n  unknowns,  where  n  is  the  number  of  interior  mesh  points. 

All  values  pertaining  to  the  boundary  mesh  points  will  be  prescribed. 
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Since  the  initial  values  for  U  for  the  interior  points  are  known  for  t  =  to,  the  solutions 
for  the  first  set  of  n  algebraic  equations  will  apply  for  t  =  to  +  p.  Successive  iterations 
will  result  in  approximate  numerical  values  of  U  for  ail  mesh  points  in  R  when  t  =  to  +  Ip. 
to  +  2p, . . .  ,to+mp  where  m  is  the  total  number  of  time  increments. 

An  explicit  procedure  for  solving  boundary  value  problems  of  the  type  containing 
Equation  (1)  is  referred  to  as  the  forward  difference  method.  Each  value  of  U  can  be 
determined  for  any  value  of  t  if  all  values  for  U  at  the  preceding  t  are  known.  Unfortunately 
the  forward  difference  method  may  be  numerically  unstable  unless  restrictive  conditions 
are  imposed  upon  a  parameter  for  selecting  p.  As  a  result,  the  maximum  value  for  p  is 
limited;  nevertheless,  for  most  cases  p  must  be  chosen  so  small  relative  to  the  space  size 
of  the  mesh  that  the  computational  labor  involved  prohibits  solutions  except  where  high¬ 
speed  digital  computers  are  available. 

Once  the  values  of  U  have  been  established  for  all  mesh  points  within  R  and  for  t  =  to+  Ip, 
to  +  2p, . .  ..to  +  mp,  then  the  rate  of  heat  flow,  q,  across  the  boundary  can  be  calculated. 

The  total  heat  removed,  Q,  can  be  determined  by  a  summation  of  q. 

m 

mp  Z  q  (2) 

I 

From  this  term,  the  total  weight  of  the  coolant  expended  is  determined  by  dividing  by 
the  enthalpy  that  is  required  to  convert  the  liquid  coolant  to  a  gaseous  state.  This  process 
can  then  be  repeated  for  a  number  of  initial  insulation  thicknesses  to  arrive  at  some 
relationship'  between  this  independent  variable  and  the  coolant  requirements. 

The  introduction  and  availability  of  high-speed  computational  machines  have  made 
practical  the  solution  of  this  problem  by  numerical  methods.  Unfortunately,  considerable 
effort  is  involved  in  formulating  the  problem  into  computer  language  and  preparing  the 
machine  program.  As  a  result,  the  time  lag  between  the  start  and  solution  of  a  problem 
might  be  significant  but  of  less  relative  importance  where  a  great  many  cases  are  to  be 
solved.  One  valuable  aid  for  expediting  problem  processing  is  an  adaptation  of  a  concise 
format  for  presenting  solution  procedures.  Section  4  contains  the  problem  statement  in  a 
recommended  form  that  can  be  understood  by  engineers  as  well  as  programmers. 


3.2  General  Assumptions 

Before  presenting  a  solution  to  the  time-dependent  heat  transfer  or  diffusion  equation, 
certain  assumptions  suggested  by  the  physical  nature  of  the  problem  should  be  reviewed. 
These  involve  the  selection  of  significant  spatial  variables,  the  manner  of  analytically 
expressing  material  properties,  and  general  definitions  pertaining  to  boundary  conditions. 


The  first  problem  can  be  resolved  by  determining  the  relative  significance  of  heat 
transfer  with  respect  to  the  coordinate  system  shown  in  Figure  10.  These  coordinates 
pass  through  an  element  of  insulation  in  an  area  of  the  vehicle  where  the  orientation  of 
isotherms  results  in  maximum  heat-flow  rates.  Since  the  thermal  conductivity  of  the 
material  is  identical  in  all  three  directions,  then  only  the  ratio  of  temperature  gradients 
in  these  directions  need  be  compared. 


Temperature  variations  with  distance  can  be  expected  to  reach  a  maximum  in  a  section 
directly  behind  the  leading  edge  of  the  wing  where  the  value  for  ^  'u^)  A  x  ”  only 
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I  Figure  10.  Element  of  Insulation 

! 

I  2  percent  of  the  heat  conducted  depthwise  through  the  structure  will  be  transmitted  in  a 

[  chordwise  direction.  In  Section  5,  you  will  note  that  the  total  weight  of  the  system  will 

k  not  be  affected  by  this  small  percentage.  The  heat  transferred  in  a  spanwise  direction  is 

I  still  smaller;  in  fact,  it  can  be  considered  as  nil. 

I  Another  section  of  the  vehicle  that  was  chosen  for  evaluating  the  merits  of  a  three- 

,  dimensional  heat-transfer  analysis  was  the  upper,  rear  portion  of  the  fuselage  where 

,  the  surface  temperatures  during  the  period  of  maximum  re-entry  heating  are  relatively 

low.  In  this  area,  the  value  of  the  ratio =  200,  which  indicates  that  the  radial 

Au  (z)  A  X 

,  heat  transfer  far  exceeds  the  circumferential  conduction;  hence,  the  latter  could  be 

neglected.  The  same  approach  can  be  applied  to  justify  neglecting  heat  that  is  transferred 
axiaily  along  the  fuselage. 

Although  we  have  shown  that  analytical  solutions  in  the  practical  sense  will  be 
sufficiently  accurate  if,  in  general,  only  one-dimensional  heat  transfer  is  assumed,  a 
numerical  approximation  of  the  diffusion  equation  in  three  spatial  variables  will  be 
derived,  since  it  provides  a  valuable  tool  for  evaluating  a  number  of  engineering 
problems  in  heat  transfer.  The  matter  of  transforming  the  derivation  to  a  one-dimensional 
equation  is  accomplished  merely  by  retaining  only  those  terms  applicable  to  the  particular 
problem  in  question. 

I  During  re-entry,  the  temperatures  of  the  external  heat  shields  may  be  related  in  some 

manner  to  mission  time.  Although  no  exact  mathematical  formulation  is  possible  for  this 
function,  it  may  best  be  approximated  by  a  series  of  points;  its  neighborhood  can  be 
j  defined  by  linear  interpolation.  This  assumption  is  not  as  wieldy  as  it  may  appear,  since 

'  digital  computer  techniques  can  readily  handle  this  procedure. 

Another  general  assumption  is  associated  with  the  mathematical  treatment  of  material 
thermal  properties;  namely,  surface  emissivity,  thermal  conductivity,  and  specific  heat. 
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The  thermal  emissivity  of  the  material  surfaces  is  assumed  to  be  independent  of 
temperature.  This  assumption  is  valid  for  two  reasons.  First,  the  emissivity  of  surfaces 
experiencing  an  appreciable  temperature  variation  (in  the  order  of  2000°F)  during  re-entry 
are  essentially  grey  bodies.  Secondly,  all  nongrey  surfaces  experience  only  a  small  change 
in  temperature,  usually  less  than  100®F. 

Thermal  conductivity  and  specific  heat,  on  the  other  hand,  are  definitely  functions  of 
temperature.  If  these  properties  are  inspected  for  a  number  of  applicable  materials,  we 
see  that  thermal  conductivity  can  be  reasonably  defined  by  a  quadratic  term  while 
specific  heat  can  be  approximated  by  a  linear  Unction.  This  approximation  results  in 
expressions  that  are  easily  operated  upon  and.  in  addition,  provides  a  convenient  means 
of  reducing  the  degree  of  the  equations  by  setting  coefficient  equal  to  zero.  For  example, 
if  the  thermal  conductivity  is  best  represented  by  a  linear  function,  a  second  degree  polynomial 
can  be  replaced  by  a  linear  function  merely  by  equatingthe  coefficient  of  the  squared  term  to 
zero.  Assumptions  concerning  the  selection  of  a  value  for  the  radiation  form  or  geometry  factor 
and  a  decision  to  om it  a ir  convection  terms  have  been  made  in  the  interest  of  generalizing 
solutions.  Since  both  of  these  terms  are  related  to  specific  locations  of  vehicle  surfaces, 
including  particular  design  considerations  would  be  necessary  for  all  cross-sectional  areas  of 
a  thermally  protected  airframe. 

Many  other  significant  assumptions  could  be  discussed;  however,  they  are  more  closely 
associated  with  particular  boundary-value  problems  and  would  not  apply  in  general.  These 
assumptions  will  be  reserved  for  Section  3.4  where  particular  boundary  value  problems 
will  be  analyzed. 

3.3  Solution  of  the  Diffusion  Equation 

The  partial  differential  equation  describing  the  time-dependent  flow  of  heat  through  the 
thermal  protection  system  is  derived  in  this  section;  it  is  then  transformed  into  an 
algebraic  expression  by  the  application  of  the  calculus  of  finite  differences.  The  final 
form  is  one  that  expresses  the  temperature  of  a  finite  element  at  time  t  +  1,  in  terms  of 
the  temperature  of  elements  at  time  t.  By  successive  application  of  this  equation  together 
with  appropriate  boundary  and  initial  conditions,  the  temperature  history  through  the 
insulation  may  be  predicted;  thereby,  the  rate  of  heat  transfer  from  the  n’^'"'  element  is 
established.  This  latter  quantity  establishes  the  cooling  weight  of  the  system  for  any 
given  amount  of  applied  insulation. 

Note  Figure  10.  The  heat  entering  face  lx  can  be  written  as 


d  u 

dQ,  =  dydz  l<  (  u  )  3-; —  dt 

Od  X 

Now  the  heat  leaving  face  2x  at  x  +  dx  can  be  determined  by  letting 


(3) 


F(x.u) 


=  k{u  ) 


dx 


(4) 


By  changing  x  to  x  +  dx,  F(x,u)  becomes  F(x  +  dx,u).  This  latter  term  may  be  expended 
by  using  a  Taylor  series,  thus 


F  (  X  +  dx,u  )  =  F  (x,  u) 


dF 

dx 


dx 


k  ( 


k 


(u) 


du 

d  X 


dx 


(5) 
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Therefore,  the  heat  leaving  face  2x  can  be  expressed  as 


'‘°x+dx 


d 

dx 


du 


k(  u  )  -  dx  1  df 
dx  d  X  ox 


This  same  approach  would  apply  to  adjacent  y  and  z  faces.  That  is 


du 


dO  =—  dxdz  k(u) 
y  dy 


dt  . 
du  ,  d 


4- dy  =  ~  dxdy  [k(u)-^  4 


dy  d  y 


k(u) 


d  u 
d  y 


dy  j  dt 


and 


du 


dQ-  =  —  dxdy  k(u  )  — —  dt  , 
‘  02 


<^024-02=-“’“'^  TT 


du  d  .,,du 

+  -  k  ( u  )  — : -  d  2  I  d  t 


(6) 

(7) 

(8) 

(9) 

(10) 


d  2  d  2  '  ’  d  2 

Since  the  net  heat  flow  into  an  element  must  be  equal  to  the  heat  stored  within  the  element, 

d0^4-  dQy  4  dOj  =  ‘‘Qx4d2  '•"‘‘Qy  4dy  +  ‘‘Q24d2  +  P Cpdx dydz-^  dt  . 

When  the  dQ  terms  are  replaced  by  the  expressions  Equations  (3),  (6),  (7),  (8),  (9),  and 
(10),  all  terms  to  the  left  in  Equation  (11)  are  cancelled  by  the  first  terms  of  the 
hrst  three  expressions  on  the  right.  The  result  is  the  following  three-dimensional 
conduction  equation: 


— -  k(u) - 

d  X 


d  X  d  y 

Equation  (11)  can  be  expanded  to 


d  d  u 

k(u) - 1-  k(u) 


dx 


dx 


d 

dz 


d  u 

d 

-  J. 

If  f  II  ) 

d  u 

'  K|  U  1 

0 

dy 

d  z 

K  V  U  / 

d  z 

d*o 
dx  * 

-  + 

— ^  k(u) 
dy 

du 

dy 

"-'-p 

d  y 

k  (u) 

d  u 

d  z 

4-  k  (u  )  ■ 

d*  u 
d 

■  -  P<=p 

d  u 


(12) 


d  u 
dt 


(13) 


Equation  (13)  can  now  be  transformed  into  an  algebraic  expression  by  substituting 
appropriate  finite  difference  relationships  for  the  partial  derivatives.  These  approxi¬ 
mations,  derived  from  interpolation  formulas,  are  listed  as  follows: 


Ot 

X  c 

-irl 

u(i4l,x,t)-u(i-l,x,t)]  , 

(14) 

du  V 

1 

2 

[u(i-t-l,x,t)  —  u(i  —  l,x,t)]  , 

(15) 

dx  ) 

4h2 

d*u 

1 

[u{i4  l,x,t)— 2u  (i.t)  +  u(i  -  l,x,t  )]  , 

dx^' 

h  ^ 

(16) 
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d  u 

1 

[u  (  i.  t-i-l  >  -  u  (  1,  f  )] 

d  t 

P 

(17) 

d  u 

1 

[u(  l+l,x,l)  -  u(  i.t  )  ] 

dx 

h 

(18) 

Although  these  expressions  have  been  written  for  the  x  direction  only,  those  for  the  y  and 
z  spatial  variables  are  similar.  In  addition,  a  square  element  was  selected  to  eliminate 
area  terms  from  the  final  equations.  The  derivation  of  these  expressions  is  included  in 
Appendix  1. 

Now  the  thermal  conductivity  of  the  material  can  be  expressed  as  a  function  of  tempera¬ 
ture  by  the  polynomial 

k(u)  =  Aq+  a,  u  +  A* 

and  the  specific  heat  by 

Cp  <**)  =  Bq  +  B,  u  . 

Taking  the  derivative  of  Equation  (19)  and  dividing  by  a  x  gives 


(19) 

(20) 


dk(u)  d  u 

'“•.r  (  A.  +  2A,U  ) 


d  It 


dx 


(21) 


Substituting  Equation  (21)  into  Equation  (13)  results  in  the  expression 


I  +  ZAj  u  (-|^)‘  -t-  Ag  +  A,  u  +  A^  0*  (-|^)  -t-  A,  + 


*  /d*  u 


2  A, 


“  dr)’-"  “o  (^1+ 


d  *  u 


d  u 


dy 


dy 


dz 


0**,"  *''<"•(-1^)  +  “I-It 


d  u 


(22) 


Replacing  tttfc  derivatives  with  the  appropriate  finite  difference  relationships  gives 

[A|  +  2AjU  (i,f)]  |[u*(i  +  l,x,t)-2u  (i+l,x,t)  u(i-l,x,f)-*-u*(i-l,x,t)]  + 

[u*(i+l,y,t)-2u{i  +  l,y,f)u(i-l,y,f)  +  u^(i-l,y,f)  ]  + 

[u*(  i+l,  z,t)-2u(i  +  l,z,f)  u  (i-l,z,t)+u^(i  -  I,  z  ,t)  ]  + 

Jz  [  Aq  +  A,  u  (  I  ,  f  ) -»-  Ag  u*  (  i,t  )j  [  u  {i  +  I ,  x,t )  — 2u  (  i,  t)  +  u  ( i  -  I,  X,  t  )  ]  + 

[u  (j  +  l,y,f)  -  2u(i,f)-Hu(i  —  (,y,t)]  +  [u(i+l,z,l)-2u(i,t)  +  u(  i-l,z, »)]]•  = 


_P_ 

P 


[  Bq  +  B,  u  (i.t  )  ]  [  u  (i.t  -I-  I)  -  u  (  i.  t)] 


(23) 
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Finally,  solving  for  u(i,t  + 1)  gives 


P  r 

u(i,»+l)  =  u(i.t)+  ,  t - T  1  A|  +  2A,u(i,t) 

20 /9h*  [Bo+  B,  u  (i.t)  ]  '■ 

[u^(i  +  l,x,f)— 2u(i-»-l,x,l)  u(i-l,x,t)-fu*(i-l,x,f)  +  u*(i+l,y,f)- 
2u(i+!,y,t)u{i-l,  y,f)  +  u*(i-l,y,»)  +  u*(i+l,7,t)-2u  (i+1,7,  t)  u(i-l,z,t)  + 
u^(  i-  I,  z,1  )]  +•  4  [Aq+  A,u(i.t)+A2U*(i.t)] 

[  u  (i  +  l,x,t)  +  u(i-l,x,t)+-u  (i  +  l,y,  f)  +  u(i-l,y,f)-f- 

u(i4l,z,  f)+o(i— l.z.f)—  6u(i,f)]  }  ^24) 


For  the  case  of  one-dimensional  heat  flow,  all  terms  associated  with  y  and  z  directions 
will  drop  out  and  Equation  (24)  reduces  to 

P 

u  (i,f  +  l)  *u(  - - - — 

20/Oh*  [Bq  +  u  (  i,t)  ] 

{  A|  +  2AgU  (i,t)  [u*(i+t,f)-2u(i+l,f)u(i-l,  t)+  u*(i-l,t)]  + 

4[Ao-i-A,  uii.tJ+AjU^ii.t)]  [u(i+l,t)-2u(i,t)+u(i-l,f)  ]}  (25) 

This  is  an  expression  for  the  temperature  of  an  element  of  insulation  of  length  h  and  of 
unit  area  at  time  t  +  1,  in  terms  of  temperature-dependent  thermal  properties  and  the 
temperature  of  the  same  elements  at  time  t. 

3.4  Boundary  Value  Problems 


A  number  of  typical  boundary  value  problems  that  were  encountered  during  the  formula¬ 
tion  of  optimization  procedures  for  thermal  protection  systems  are  presented  in  this 
section.  Each  problem  is  associated  with  some  definite  combinations  of  thermal  insula¬ 
tion  and  cooling  system  or  thermal  insulation  without  a  cooling  system,  which  are 
representative  of  schemes  evolving  from  current  technological  developments.  The 
following  types  will  be  analyzed  in  the  following  sections:  insulated  and  cooled 
structure,  insulated  and  cooled  compartments,  and  insulated  structure. 


3.4.1  Insulated  and  Cooled  Structure 

A  cross  section  through  the  wall  of  an  insulated  and  cooled  structure  is  shown  in 
Figure  11.  Let  s  represent  the  external  surface  of  the  radiation  shields,  o  the  interface 
between  the  heat  shield  and  thermal  insulation,  r  the  inboard  surface  of  the  insulation. 
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c  the  cooled  structure  and  g  the  compartment  atmosphere  and  insulation  interface. 


HEAT  SHIELD 


, - ^  PRIMARY  , - ^  STRUCTURE  , - ^ 

>  INSULATION 

COMPARTMENT  ATMOSPHERE 


Figure  11.  Insulated  and  Cooled  Structure 


The  conditions  of  the  boundary  value  problem  are  stated  as  follows: 
Considering  the  outer  layer  of  insulation  (o  to  r)  gives 


Ut 

(  X  ,  t  ) 

=  G  (U  ) 

0  <  X  <  r  ,  t  >  0  , 

(26) 

U, 

(  X  ,t  ) 

=  [u(s,t)-u(o,t)l 

k(u)  '■  ■' 

X  =  0,  t  >  0 

(27) 

Ux 

(X  .  t) 

=  [u(r,t  )-  u  (c)  ] 

k(u)  *■  ■' 

M 

II 

V 

o 

(28) 

and 

U(k,o)  =  constant  0<x<r,t=0- 

A  similar  set  of  conditions  applies  to  insulation  contained  by  the  cooled  structural 
shell.' 


U(  ( X, t )  =  G  I  u ) 

U„  (X  ,f  )  =  f  u  (r,  t)  -  u(c)] 

“  k(u)  J 

[u(g.t)-  u  (V)] 

U  (x,o)=  con  St  ant 

The  problem  statement  can  now  be  completed 
affect  a  solution  thus 

U  (s  ,t)  =  F(t  ) 

U  (c,t)  =  constant 

22 


c  <  X  <  g  ,  t  >0  , 

(30) 

X  =  c  ,  t  >  0  , 

(31) 

X  =  g  ,  t  >  0  , 

(32) 

c<x<g,t  =  0- 

(33) 

stating  the  remaining  relationships  that 

X  =  s,  t  >  0  , 

(34) 

o 

Al 

u 

II 

X 

(35) 
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U(v,f)=constant  x=v,t>0  (36) 

The  preceding  conditions  can  now  be  transformed  into  algebraic  expressions  by  replacing 
the  derivatives  with  the  corresponding  finite  difference  relationships  listed  in  Section  3.3, 
Equations  (14)  through  (18). 

The  transformation  of  the  identical  Equations  (26)  and  (30)  was  discussed  in  Section  3.3. 
As  a  result,  Equation  (25)  was  derived.  Similar  operations  on  Equation  (27)  will  now  be 
performed; 


k(u) 


d  u 


=  C, 


[u(s,t)  -u(o,f)] 


(37) 


where  Cjj  is  the  conductance  of  the  radiation  heat  shield.  Referring  to  Figure  12, 


Cb  =  +  C„  -t-  C,  (3g) 

states  that  the  overall  conductance  of  a  honeycomb  sandwich  panel  is  the  sum  of  the 
conductance  through  the  metal  core,  Cj.  and  the  conductances  of  the  air  space  C  +  C  . 
Now  ’’ 


Cd 


“d  S,,  ^  ka  (V  )So 

^  'O'  I 


Cr  =  c/iu)  So 


(39) 


But  Cb  can  be  approximated  by  a  quadratic 

=  Do  -H  D,  u  + 


(40) 


completing  the  transformation 

"T  (  +  (o+l,»  +  l)]-f  [u*  (  o.t  +  I )  + 

2u(o,f+l)  u(o4-l,t+l)  +  u*(o  +  l,t+l  )]}[u(o,t  +  l)-u(o  +  l,f+l)]  = 
Do  +  [u(s,f  +  l)  +  u(o,t+l)]+-^  [u*(s,l  +  l)  + 

2u{s,t  +  l)u(o,t+l)  +  u*(o,t  +  l)][u  (s,t+l)— u(o,t4l)^ 


(41) 
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Expanding,  collecting,  and  transposing  terms.  Equation  (41)  becomes 

(  Aj+Oj  h)  u  ’  (o  ,  f  +  I  )  +  [2A,  +Aju  (  0  +  I,  t  +  I  )  +  2hD,  +  hDjU  (s.  t  +  I  )]u*(o  ,t4  l)  + 

[4hDp-  h  Da  u*  (s  ,  t  +  I  )  -  A^a^io  +  I,  t  +  t)+  4  Aq]  u  (o  ,  t  4  I  )  -  [  AAqu  (o  4  I ,  f  +  I )  4 

2A,u*(o4l,t4l)4AjU*(o4l,t4l)4  4hDou{s,t4l)42hD,  u*(s,t4  I  )  4- 

hD2U*(s,t4l  )j  =  O*  ^42) 

The  term  u(s,t  +  1)  =  F(t)  is  given  and  u(o  +  l,t  +  l)  can  be  determined  from  Equation  (25); 
therefore,  u(o,t  +  l)  can  be  found  by  solving  Equation  (42). 

The  transformation  of  Equation  (28)  is  accomplished  in  a  similar  fashion 

k(u)-^l  =C„  4C,  [uir.  t)  -  u  (c)]  ,43. 


where 


— — -  =  — —  (A  4  A„  u_  +A--U*) 

j  /  '  00  01  fn  02  m 

-^o  ^o 


where  is  the  thermal  conductivity  of  air. 


In  Equation  (44), 


u{r,t4l)4  u(c) 


Cr  =  ^  >^0  { 


Expanding  Equation  (46) 


[459.  69  +u(r,t4l)]^-[459.69  4u(c)] 


u  ( r-,  f  4  I )  -  u  (  c  ) 


(T  F, 


Cf  =  - - — 5 - r  (  r  ,  t  4  I)  4  1838  eu®  (  r  ,  t  4  I  )  4 

u(r,t4l)-u(cl  •• 

126  7.9  X  10®  u*  (r.t  4  I  )  4  388  5.6  X  I O®  u  (r.t  4  I  )  - 

u'*  (c)  -I838.8  u®{c)  -  1267.9  XIO®  u*(c)  -  3885  6  X  I0®u(  c)] 


Substituting  Equations  (44),  (45),  and  (47)  into  Equation  (43)  together  with  the  finite 
difference  expression  for  the  partial  derivative  yields  the  following: 
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u(r,t4-l  )‘l-u(c) 


T  r  u  (  r,  t  -f  I  )  +  u  (c)  it  ■, 

]  +  [ - i - ]  } 


ot 


u  (  r  ,t  4  I )-«- u  (c)  1*  \  r  u(r,t  4  I) -u  (c) 


]  }( 


-]  (r.  »+')<• 


18  38.8  <r  u*  (  r ,  I  4  I  )  4  1267.  9  X  10*  <r  F,  F,  u*  (r .  t  4  I  )  4 

3885.6  X  10* <r  F^  F^  u  (r.t  4-1  »-  <r  F,  F^  (c)  -  1838  8  aF^  F^  u*  (c)  - 

1267.9  X  to*  trF^  F„  o*(c)  —  388S.6  X  10*  «r  F,  F,  u(c) 

Equation  (48)  reduces  to 

4h«rF^F^u^  (r,t4  1)4*  +  7355.2  h  «r  F^  F„  -  A,  )u*  ( r,  t  4  I  )  4 


(48) 


„(c)  +5071.6  X  IO*h«r  F*  Fo  -2  A,  -A,  u(c)]u*(r,t4l)4 
^  4^Ago  ^  +i554  2xlO*ho-F,  F^  -4Ao-  u*(c )]  u  (r,l  +  I )  + 


4Aou(c)  +  2A,u*(e)  •+Aju*(c) - °**  u(e) - u*(c)  - 

A  O 


2h  A, 


■**  4---  u*(c  )  -  4h<r  F^  F„  u*(c)  -  73  55,2  h  cr  F^  Fg  u»(c)  - 

-^o 

5071.  6  X  lO’  ho- F^  Fg  u*(c»  -  15542  X  10*  h  (T  F^  Fg  u(c  )  ]  s  O 

Again,  this  equation  can  be  solved  for  u(r,t  + 1),  since  all  other  terms  are  known  or  can 
be  determined  by  previously  established  relationships; 

We  assumed  by  Equation  (35)  that  the  circulating  coolant  media  maintains  a  constant 
structural  temperature;  hence,  the  heat  capacity  of  these  elements  does  not  affect  the 
temperature  gradient.  Furthermore,  the  thermal  conductivity  of  common  materials  used 
in  this  area  are  very  high  as  compared  with  those  of  thermal  insulations  and,  when 
coupled  with  typically  thin  sections,  the  temperature  drop  across  the  member  is 
negligible. 

Continuing  with  the  solution,  the  temperature  of  internal  elements  on  the  Inboard  side 
of  the  structure  can  be  determined  by  Equation  (30).  Again,  the  transformed  relationship 
is  identical  to  Equation  (25). 

Now  Equation  (31)  states  that  the  heat  transferred  from  r  to  c  is  equal  to  the  heat 
absorbed  by  the  flrst  element  of  insulation  adjacent  to  c.  A  solution  in  this  case  was 
previously  obtained  (Equation  (49)),  since  the  temperature  of  c  is  constant;  therefore, 
u(c)  is  independent  of  conditions  beyond  that  station. 
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The  relationship  remaining  to  be  solved  is  Equation  (32)  thus 

k(u)  -|^  [  =  [uig.f )  -  u(v)]  (50) 

»=  9 

Again,  the  partial  derivative  is  transformed  to  the  difference  relationship 
"h"  {*0+  +  4^  [u*(g-l,t+l)  + 

2u  (g-l,  t  +  I  )  u(  g,t  +  I  )  -Hu*  {g.t  i-l  )]  }  [  u  ( g-l ,  t I )  -  u  (g  ,  f  +1  )  ]  = 

%  [u  (g.t-HI  )  -u(  V)] 

Expanding,  grouping  terms,  and  transposing.  Equation  (51)  becomes 

Aju*(g,f +  1  )  +  [2A,  -  Aju  (g-  I  ,  t  -H)]  u*  (g,t+l  )  ^[aAo-Aj  u*(g-l,  ti- 1) +4hcju(g,t+l)- 


[4Aou(g-l.  t  +()-»-  2A,u*(g-l,t+l  J  +  A,u*  (g-l  ,  t  +1  )  +  4hC^  u(w)  ]  = 


(52) 


The  solution  of  Equation  (52)  for  u(g,t  + 1)  completes  the  process  of  determining  the 
temperature  distribution  through  the  thermal  protection  system  at  time  t  +  1.  A  repeti¬ 
tion  of  the  process  results  in  the  temperature  distribution  at  time  t  +2.  By  this  technique, 
the  temperature  distribution  for  all  t's  may  be  established. 


An  internal  cooling  system  maintains  u(v)  constant  by  circulating  the  compartment 
atmosphere  across  the  face  of  the  insulation.  A  similar  condition  existed  at  c  where  a 
fluid  is  circulated  through  passages  that  are  attached  to  the  structure.  Determining  the 
weights  of  coolant  in  each  case  requires  that  the  rates  of  heat  removal  at  c  and  v  must 
first  be  calculated. 


At  the  structure,  c, 

q(c,f  +  l)  =  Cg+C,  [u(r,t-tl)-u(c)] 


(53) 


since  u(c)  is  a  constant.  Equation  (53)  simply  states  that  the  rate  of  heat  transfer  to  the 
coolant  at  time  t  +  1  equals  the  rate  of  heat  transfer  across  the  air  gap  at  time  t  +  1.  The 
total  weight  of  coolant  expended  during  re-entry  is  then 


W, 


fnp 


'f9 


Z  q(c,t). 

I 


(54) 


Novs^the  rate  of  heat  transfer  to  the  compartment 

q|  V,  t-HI  )  =  Ch  [u  (  g,  t  +  I  )  -u  (v)] 


and  the  weight  of  expended  coolant 

m  p 


Z  q  (  V,  t)  . 


"fg 


(55) 


26 


ASD-TDR-62-625 


The  amount  of  insulation  in  each  case,  which  results  in  a  minimum  combined 
insulation  and  coolant  weight,  can  be  determined  by  repeating  the  solution  for  several 
insulation  thicknesses;  then  the  combined  weights  can  be  compared  with  the  Insulation 
thickness. 

3.4.2  Insulated  and  Cooled  Compartment 

A  typical  section  of  this  construction,  consisting  of  an  insulated  radiation  shield 
exposed  to  the  airstream  and  an  insulated  cooling  system  attached  to  the  compartment 
wall,  is  shown  in  Figure  13.  The  various  surfaces  and  interfaces  are  identified  as 
follows: 

a.  At  X  =  s,  the  exterior  surface  of  the  radiation  shield  is  supported  from  a  stiffened 
panel,  f,  by  refractory  metal  members. 

b.  At  X  =  o.  a  retaining  sheet  confines  thermal  insulation  in  an  area  that  is  adjacent 
to  the  cooling  system,  r.  This  assembly  together  with  another  layer  of  insulation  is 
supported  from  the  compartment  wall,  c. 

c.  The  compartment  atmosphere  is  located  at  x  =  v. 

The  boundary  value  problems  associated  with  this  case  can  now  be  specified. 
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Figure  13.  Hot  Structure  -  Insulated  and  Cooled  Compartment 


Conditions  tor  the  insulated  radiation  shield  are  as  follows: 


U,  (X,  t ) 

=  G  (u) 

s<x<f,t>0  , 

(56) 

U  (s.t  ) 

=  F,  ( t  ) 

X  =  a  ,  t  2  O  , 

(57) 

U^(x,f)  = 

;^{(Co  +  (V)[u«f.t)-u(o,t)] - 
Wu)  ^  ■* 

C,  [u(s,l)-u  <».t  )]  } 

X  =  f  ,  t  >  O  , 

(58) 
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U(x,o)scon3tant  s<x<f,t  =  0. 


For  the  insulation  between  x  =  o  and  x  =  r. 


U,  (x,t)  r  G  (u) 

0  <  X  <  r,  1  >0  , 

(60) 

X  =  0  ,  t  >  0  , 

(61) 

X  =  r  ,  t  >  0  , 

(62) 

LI(x,o)  =  con  St  ant 

o<  x<r,  t  =  o. 

(63) 

Finally,  for  the  remaining  section  of  insulation  between  x 

=  r  and  x  =  c 

U,  (x.t)  =  G(u) 

r  <  X  <  c  ,  t  >0  , 

(64) 

X  =  r,  t  >  o  , 

(65) 

U,(x,t)  =  — ^  f*l.  (u  (c  ,  t  )-u(w))  1 

X  L  s  B  \  '  J 

x  =  c  ,  t  >  o  , 

(66) 

U  (  X  ,0 )  =  c  on  stonf 

r  <  X  <  c,  t  =  0 

(67) 

Relationships,  in  addition  to  those  just  listed,  which  will  be  required  for  a  solution, 
as  follows: 

are 

U  (r  ,t  )  =  (  t  ) 

K 

II 

IV 

o 

(68) 

►'fg  (r)  =F,  (t  ) 

X  =  r ,  t  >  o  , 

(69) 

U  ( V  ,  t  )  =  con  stonf 

X 

II 

< 

IV 

O 

(70) 

Following  the  procedure  applied  in  Section  3.4.1,  the  partial  derivatives  of  the 
boundary  condition  are  transformed  to  finite  difference  expressions.  Equations  (56),  (60), 
and  (64)  reduce,  as  before,  to  Equation  (25). 


Since  the  temperatures  of  x  =  s  for  all  t’s  are  given,  and  the  temperature  of  internal 
elements  of  insulation  for  the  first-time  step,  t  + 1,  can  be  determined  from  the  solution 
of  Equations  (56)  and  (60),  then  the  next  operation  should  produce  the  temperatures  at 
x  =  f  and  x  =  o  for  this  same  time  step.  Now  the  rate  of  heat  transfer  from  the  insulation 
element  at  x  =  f  must  be  equal  to  the  rate  of  heat  transfer  to  the  insulation  element  at  x  =  o 
if  the  material  in  the  region  f  £  x  <  o  has  zero  heat  capacity.  This  being  the  case  for  air, 
the  two  applicable  Equations  (58)  and  (61)  at  this  location  can  be  written  as  follows: 

k(u)-j^|  +  C,  [u  {  s.t  )  -  u  (  f.f  )]  =  (Cg  +  C,  )  [u(f,»)  -  u(o  ,n], 

X=  f 
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The  third  coefficient  is  equated  to  the  thermal  conductance  of  the  air  space  established 
by  Equation  (44),  in  which  Um  is  now 


u(f.t)  +■  u  (o,t) 


(81) 


Therefore, 


Co  =  C, 


[4Ao„+2Ao,  uif.f)  +2Ao,u(o,t)  +  Ao^u^if.t)  + 


2Ag2u(f,t)  u  (o,t  )  +  Aq2  u  io.f)- 


(82) 


Considering  the  relationship  for  the  equivalent  thermal  conductance  for  air-gap  radia¬ 
tion,  an  expression  similar  to  Equation  (47)  is  obtained.  Rewriting  Equation  (47)  and 
inserting  the  proper  subscripts  gives 


Cr  = 


c^  I  f  r  4 

-  =  -  {  O-  ^  u‘’(f.t)  + 

u  (f.  t  +  1)  -  u  (0.  f  +  I )  o(f  ,1  +  1)- u  (o,f+l)  *•  «  Oi¬ 

ls  38.  8u’(f  ,t  )+  1267.9X10*  u*  (f,f)  +  388  5.6  X 10*  u  (f  ,1)  - 

u'*{o,  t  )  -  1838.  8  u’  (o.f  )  -  12  6  7.  9  X  lo’  u*  ( o,  t  )  -  3885.  6  X  lO*  u  (o,f  )]} 
Equation  (83)  shows  that  C4  is  equal  to  the  term  enclosed  by  the  backets. 


(83) 


The  final  coefficient  is  derived  from  the  thermal  conductance  term  for  the  element  of 
insulation  at  x  =  o  by  a  manner  similar  to  that  utilized  in  establishing  Equations  (77) 
through  (80)  where  is  now 


u  (o,t  )  +u  (0  +  1 ,  t+  I ) 


(84) 


Now 


I  r  A| 


u(o+  I,  t  +  1  )+-^  u*(o,t  )  + 


U  (o,t)  U  (0+1, t  +  1)  +  ^  (o  +  l,t+l)j 


(85) 


Recalling  the  simultaneous  Equations  (71)  and  (72),  the  partial  derivatives  may  be 
replaced  by  finite  terms  and  substituting  the  coefficients  for  the  conductance  and  thermal 
conductivity  results  in  the  following: 


[u(f- I,  f +1)  -  u  (  f ,  t +1)]  -l-C,  [u(s,t +1)  -u(f,t +  1)  ]  = 

C,[u  (f,  t+1)  -  u(o,  t-l- I  )  ]  +  C4  , 

C,  [u(f,t  +-l)-u(o,t+l)]+-C^  =  [  u{o,  t  +1)  -u(o+  l.f  +  I  )  ] 


(86) 

(87) 
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Rearranging  Equations  (85)  and  (86)  by  transposing  terms  gives 

(c,  +C,-»-  -^)  [u(f,i  +  l)]  -C,u(o.t+l)  = 

u(  f  -  I,  t  +  I  )  +  C, 

c,  u  ( f.t  +  n-(c3+  -^)  u(o,t  +n  =-C4  - 

The  temperature  u(f,t  +  l)  and  u(o,t  +  l)  can  now  be  determined  from  Equations  (88)  and 
(89). 

The  next  step  involves  determining  the  temperature  of  internal  elements  of  insulation 
for  r  <x  <c  at  time  t  +  1.  Following  the  procedure  employed  for  each  of  the  preceding 
insulation  thicknesses,  these  temperatures  are  determined  from  the  expanded  form  of 
Equation  (64), 

When  conditions  at  x  =  r  are  explored,  several  possibilities  may  be  encountered.  These 
involve  a  knowledge  of  the  rates  of  heat  transfer  about  and  at  the  junction  r  as  specified  by 
boundary  conditions  in  Equations  (62)  and  (65).  Recalling  that  x  =  r  is  the  location  of  the 
cooling  system  implies  that  a  heat  sink  exists  at  this  point  and  the  rate  of  heat  dissipated 
by  the  coolant  is  denoted  by  qj.  Furthermore,  the  materials  used  in  the  construction  of 
the  cooling  system  can  store  heat  at  a  rate  equal  to  qg.  Now  the  heat  involved  in  both 
mechanisms  must  be  the  balance  of  that  transferred  to  the  junction  through  the  Insulation 
at  the  rate  qj  and  the  rate  of  heat  transfer  from  the  junction  qo  so  that  the  expressions 
for  each  can  be  written  and  acted  upon  accordingly,  thus 

Oj  *  Mu)  -|^|  *  IT  {*0+  [u{r-l.t+l)  +  uir.t+l)]  + 

-r 

[u*(r  -t.f  +  I  )  4-2u(r-l,t  +  l)  u  (r,t +  1)  +  u*  ir,t  +  I )]  } 

[u(r-l,t  +  l)  -o(r,t+l)] 
or 

q,  =  +  I)  -4AoU  ir,l+l)  +  2A,u*(r-l,l+l)-2A,u*{r,f4-l)  + 

Aju’ir  -l,t  +l)+  AjU*  (r  -( ,f  +  I )  u(r ,  t4- 1)  -  A^u  (r  -  I,  t  +  I)  u*(r,  t+l)-AjU*(r,  1 ^9^^ 

In  a  similar  manner 


u(  »,H-I  )  -  ,  (88) 

O5 

-jp-  u  ( o  + 1 ,  f  ) .  (89) 


=  J_  (4A„u  (r.t  +  l)-4A„u  (r  +  I,  t  +  I  )  +  2  A,  u*(  r ,  t  +1  ) - 
4h  °  ' 

X  — +r 

2  A,  u*(r  +  I,  f  +  1  )  +  Aju’ir.f-I-I  )  +  A^u*  (f.t  +  I  )  u  (  r  +  I,  f  +  I  )  - 
AjU(r,  t+l)  u*(r +1,1  +  1  )  —  AjU*(f  +  l,t+l)  ^ 


(92) 
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The  heat  stored 


du 


at 


ft-l 


■(  p,  S/'®'  Pz\^z 

/),  Cp  ASs)!  u<r.t  +  l)-u(r,t)] 


(93) 


The  subscripts  of  the  p.  Cp,^,  and  S  terms  refer  to  particular  components  of  the  cooling 
system. 

Now  the  three  conditions  that  apply  at  the  junction  x  =  r  are  as  follows: 

a.  qj  =  +  q^.  then  q^  =  0  (94) 

b.  ^li  ■  =  '’d 

c.  <1;  qg  implies  that  u(r,t  +  l),  although  given,  is  incorrect  since  it  leads  to  an 

''impossibility  equation.”  The  solution  is  to  determine  a  new  u(r,t  +  l)  that  will  satisfy  the 
equality  qi  =  q^  +  qs-  One  must  understand  that  cy  =  0  for  this  condition. 

Combining  Equations  (91).  (92).  and  (93)  and  solving  for  u(r,t*l)  gives 

u»(r,  t  +  l)-(-  [u(r-l.t +l)+u(r+l,i  +1)]}  u*(r,f4-l)  + 


{[t  ip,^p,A  +  + 


4  h 


^  u*(r-l,f  +  l)- 


u®  (r  +  l,  14-1)}  u  (r,f +  1)  =  u(r-l,  t  4-  1)4--^  u*  (  r -I ,  t  4- I  )  4- 

4  h  ■*  h  2  h 

u’(r-  I  ,f  4  I )  4-  u  (r4-  I,  f  4-n  4-  “  (r  4-  I,  f  4  I)  4- 

^  h  n  Z  n 

TT  + [f(^.Cp  S,  4-p^C  4  i,  S3  )]u(r.t) 


(96) 


The  only  temperature  remaining  as  an  unknown  is  the  temperature  at  x  =  c.  The  boundary 
condition  for  this  location  is  defined  by  Equation  (66)  or  in  slightly  different  terminology 


k(u) 


ax 


■  Pc  ^P- 


X  =c 


Cji  [u(c,t4l)  -u{v)] 


(97) 


14  I 


Transforming  the  derivatives  to  finite  difference  terms  and  substituting  in  the  expression 
for 


k  (u)  =  Ap+  A,  u„4  AjU^* 


(98) 
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where 


u(c-l,t+l  )  +  u  (c,t) 


(99) 


results  in  the  following  relationship 

A,  .  .  A 


I  r  A|  A|  A? 

—  A„+— ^  u  (c-l,f-H  )+  — f-  u  (c.t  )  +  — ^  u*  (c-l,  t  +  l  )  + 

h  L  o  2  2  4 

-|^u(c-l,t+l)u(c,t)  +  u^'ic.t  )]  [u(c  -  I  ,t  +  l)-u  (  c  ,1  +1)  ] 


^  Pc  “^4  r  1  r  1 

-  ^ -  [u(c.f+l)-u(c.»)J  +  Ch[u(c,f-H)-u(v)J 


(100) 


To  similify  the  algebra,  let 


C.  =  [Ao-(-  -§!-  u(c-l,t  +  l)  +  u(c,f)  +  -^  u%-l,f  +  l)- 


'I  h  L-'o^  2 
A 


uCc-l,*+n  u(c,l)  H - u*(c,l)l 

4  J 


and 


Then 


^Pc  ‘'Pc  '^c 


u  {  c  ,  f  +  I 


— T  Tt;  “  [c,  u(c-l,  t+n  +  C^u  (c,t)  4-C|,u(v)  1 

+  C,  +  Ct,  '•  J 


(101) 


(102) 


(103) 


For  this  case,  the  weight  of  the  coolant  expended  during  re-entry  is  determined  by 
summing  up  the  as  determined  by  Equation  (95). 

m  q(r,f) 

=  ™>P  Z 


hfg^ 

The  heat -transfer  rate  to- the  compartment  is 

q{v,»)  =  C^[u(c,r)— u(v)j 
and  the  weight  of  the  coolant  required  at  this  location  is 

m  q  (v, t  ) 

*•  ■  ""  ?  ■ 


(104) 


(105) 


(106) 


The  combination  of  insulation  and  coolant,  which  results  in  a  minimum  system,  weight, 
is  determined  as  described  in  Section  3.4.1. 
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3.4.3  Insulated  Structure 

In  some  areas  of  re-entry  glide  vehicles  such  as  wing  surfaces  where  aerodynamic 
heating  is  much  more  intense  on  the  lower  surface  than  on  the  upper  and  the  hotter 
sections  internally  view  the  colder  ones,  heat  may  be  dissipated  through  the  structure 
by  radiation  and  gas  conduction.  This  principle  is  illustrated  in  Figure  14. 


HEAT  SHIELD 


j 

INSULATION 

_J 

STRUCTURE 

HEAT  SHIELD 

- 

Figure  14.  Hot  Structure  -  Insulated 


At  X  =  s,  the  bottom  surface  of  the  wing,  exposed  to  the  air  stream,  is  insulated  to 
limit  the  temperature  of  the  internal,  load-carrying  structure  at  a  level  where  the 
stiffness  is  not  significantly  affected.  Heat  absorbed  by  the  surface  at  x  =  s  is  conducted 
through  the  insulation  and  the  heat  shield  support  to  the  surface  at  x  =  f.  From  this 
surface,  it  is  transferred  by  radiation  and  gas  conduction  to  the  upper  surface  at  x  =  o, 
which  is  at  a  lower  temperature  because  of  the  small  heat  input  from  the  boundary 
layer.  The  problem  in  this  case  involves  the  application  of  sufficient  Insulation  to  obtain 
the  desired  temperature  at  x  =  f. 

The  conditions  for  the  boundary  value  problem  may  be  stated  as  follows: 

In  the  region  of  the  insulation 

U,(x,i)=G(u)  s<x<f,  t>o, 

( X ,  I  )  =  — ^ *  =  ».  *>0, 

k(u)  9  ^ 

where  the  heat  generated  at  the  lower  surface 

qg  =  F^(t  I  X  =  s,  f  >0  , 

(K,  f  )  =  -!— •  {  (  Ca  +  Cf  )  [u  (f  ,t  +  n-u  (o.t  +  1)1- 
^  k(u)  ■* 

C5[u(s.t+n  -  U  (f.  t  +  I  )]  }  X  =  f  ,  t>  o  , 

U  (  X.  o)  =  constant  s  <  x<  o,  t  =  o.  (m) 


(109) 

(110) 


(107) 

(108) 
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The  equation  that  applies  at  x  =  o  is 

(Co  +  Cr)[u(f,t  +  l)-u(o. !  +  !)]  =  -Tg  X  =  o,  t>o 

Here  qg  is  the  heat  generated  at  the  upper  surface  where 
Tg  =  ^2  ( i  )  X  =  o  ,  t  >  o  • 


(112) 


(113) 


After  first  assuming  some  insulation  thickness,  the  temperatures  within  the  insulation 
for  time  t  +  1  are  calculated  after  first  transforming  Equation  (107)  to  the  form  in 
Equation  (25). 

Now  the  temperature  at  x  =  s  for  time  t  +  1  can  be  determined  from  Equation  (108) 
where 


k  (u) 


dx 


=  ‘’g  - 


(114) 


The  term  q^.  is  that  heat  radiated  to  space  from  the  surface  of  the  heat  shield. 
Neglecting  the  temperature  of  heat  sink  in  space  gives 

xo-F^Fg  [u^(s.»-H)  +  l838.8u*(s.f  +  l)+l26  7.9XI0’u®(s,f-H)-l- 
3885  €  X  10®  u  (S  ,  »  +1  )  +  44  65.  4  X  lO^  ] 


and 


where 


k(u)  =  Ag  +  A,  -f  A^ 


_  u(s.t-FI)+u(8+l.t4-l) 


(116) 


(117) 


Equation  (108)  can  be  written  as 


[u(s.t  +  l)  +  u(s  +  l,t  +  l)]+-^  [u(s,t+l)  +  u(s  +  l,t  +  l)]*  } 

[  u(s,  t  +  I  )  -  u  (s  i- l,f  +  n]  =  -  <7-  Fg  Fjj  [u^  (s.t-l- I  )+  1838.8  u’ (s  ,t  +  I)  + 

126  7.  9  XI  o’  u*(s,t  +  l)  +  3885.6X10®  u(5.t+l)  +  44  65.4  X  10^] 


(118) 
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4.0  DIGITAL  COMPUTER  SOLUTIONS 


4.1  General 

Before  problems  of  the  type  outlined  in  Section  3.4  can  be  solved  by  digital  computer 
techniques,  a  concise  mathematical  statement  must  be  formulated  outlining  all  program 
requirements.  The  purpose  of  this  statement  is  to  give  the  programmer  an  insight  into  the 
type  and  scope  of  the  problem  to  be  solved  by  describing  the  algebraic  relationships 
approximating  the  set  of  differential  equations. 

The  format  for  this  section  includes  the  construction  and  usefulness  of  networks, 
required  program  inputs  and  outputs,  the  solution  procedure,  and  a  program  flow  chart. 
The  type  analyzed  in  Section  3.4.2  will  be  used  to  illustrate  the  manner  of  presenting  the 
problem. 

4.2  Mesh  Networks 

Once  a  configuration  for  the  thermal  protection  system  has  been  established,  an  array 
of  regularly  spaced  straight  lines  that  are  mutually  perpendicular  to  the  heat-flow  path 
are  constructed  within  the  boundaries  starting  from  some  arbitrary  but  fixed  point  u(Xo, 
to),  and  positive  number  h,  the  mesh  size  for  the  variable  x.  In  addition,  a  time  mesh 
size  p  is  chosen  and  the  mesh  consisting  of  another  set  of  parallel  lines  constructed 
normal  to  the  t  coordinate.  The  method  of  selecting  an  hand  pis  discussed  in  Section  4.4. 

The  network  for  the  problem  in  Section  3.4.2  consists  of  three  mesh  as  shown  in 
Figure  15;  one  for  each  section  of  thermal  insulation  at  s  <x<f,o<x<r,  and  r  <  x  <  c 
for  t  >  0.  Interior  mesh  points  are  those  lying  within  the  boundary  mesh  points  at  t  =  o, 
and  X  =  s.f,  o,r,  and  c. 

Now  the  temperatures  at  all  interior  mesh  points  can  be  found  by  solving  the  appropriate 
difference  equation  mn  times,  where  mn  is  merely  the  total  number  of  interior  mesh  points. 
Boundary  mesh  point  temperatures  are  either  given  or  can  be  calculated  from  the  relation¬ 
ships  which  apply  at  a  particular  interface. 

• 

4.3  Program  Input 

This  particular  portion  of  the  problem  statement  includes  all  known  data  that  are 
required  to  obtain  a  solution.  Referring  again  to  the  case  in  Section  3.4.2,  this  informa¬ 
tion  can  be  listed  as  given  in  paragraphs  4.3.1  through  4.3.11.  If  sufficient  data  are 
involved,  presentation  in  tabular  form  might  prove  advantageous. 

4.3.1  Surface  Temperatures 

The  temperature  of  the  surface  is  given  as  a  function  of  time  by  listing  a  sufficient 
number  of  temperatures,  u(s),  and  related  times,  t,  toadequately  represent  this  function. 
Since  an  exact,  concise  mathematical  formulation  is  usually  not  possible,  interpolation  can 
be  employed  to  approximate  values  for  u(s,t),  which  lie  between  given  points.  Example; 

Time  increment,  A  t  =  1  minute 

Total  flight  time,  tj-  =  200  minutes 

Surface  temperatures,  u(s)  =  SO.O'F,  85.0”F,  90°F  .... 
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U„. 


u, 


03  3 


“'034 


0 

^110 

UflO 

^410 

^*20 

^320 

U420 

^030 

U,40 

“240 

u,„ 

Uh, 

u,., 

U4I. 

U,2, 

U22. 

U,2I 

U42I 

Uo«i 

U|4I 

^31* 

U,*2 

U222 

*^32  2 

^42  2 

U„3 

^*13 

^313 

U4.3 

L 

U,23 

U223 

U3.4 

IJ4.4 

p 

_ , 

T 

h 

u 

iim-i 

u 

21m -1 

u 

5im-  1 

u 

41  m-i 

u 

II  m 

U 

3im 

U 

4ini 

1 

J 

12  m 

U 

22m 

‘J32m 

u 

42m 

u 

06m 

u 

14m 

u 

24m 

Figure  15.  Mesh  Network 


4,3,2  Initial  Temperatures 


Since  three  mesh  networks  are  necessary  to  depict  the  approximate  solution  for  the 
case  in  Section  3,4,2,  a  general  temperature  term  can  be  introduced  denoting  the  space 
and  time  position  associated  with  a  given  or  calculated  temperature.  Therefore,  u(i,j,m) 
identifies  a  particular  interior  mesh  jx)int  where  i  refers  to  the  spatial  mesh  column  in 
any  mesh,  j  is  a  particular  mesh,  and  m  the  time  position  or  row  in  the  mesh. 

The  initial  temperatures  can  now  be  written  as  u(i,  j,o)  and  are  given  for  all  i.j, 

4,3.3  Spatial  Mesh  Size 


Sufficient  accuracy  will  usually  be  obtained  if  the  spatial  mesh  size,  h,  is  made  small. 
Therefore,  the  number  of  interior  mesh  points  and  the  number  of  equations  to  be  solved 
will  generally  be  large.  Some  compromise  is  involved  in  selecting  a  value  for  h  that  will 
permit  satisfactory  accuracy  with  a  reasonable  number  of  interior  mesh  points. 

4.3.4  Temperature  at  x  =  r 


The  temperatures  of  the  spatial  mesh  column  x  =  r  are  equal  to  the  saturation 
temperature  of  the  coolant  used  as  the  expendable  heat  sink  and  is  a  function  of  ambient 
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pressure,  flight  altitude,  or  the  pressure  within  the  system.  Example: 


For  0  <  t  <  150; 

Time  increment,  A  t  =  150  minutes 
Total  time,  tr  =  150  minutes 
Temperature,  u(r,t)  =  80®F 
For  150  <  t  <  200: 

Time  increment,  A  t  =  5  minutes 

Total  time,  tj.  =  50  minutes 

Temperature,  u(r,t)  =  80.0"F,  85“F,  90“F  .  ,  . 

4.3.5  Densities 

The  density  of  each  material  used  in  the  construction  of  the  thermal  protection  system 
is  as  follows: 

PyP2.Py  ... 

4.3.6  Thermal  Conductivities 

The  thermal  conductivity  of  each  material  used  in  the  construction  of  the  thermal 
protection  system  as  a  function  of  temperature  is  also  required.  Example: 

J  =  1: 

u  =  100”F,  200”F.  300'F  ,  ,  . 
k:(u)  =  1.0,  1.1,  1.2  ..  . 

4.3.7  Specific  Heat 

The  specific  heat  of  each  material  used  in  the  construction  of  the  thermal  protection 
system  as  a  function  of  temperature  is  also  required.  Example: 

J  =  1: 

u  =  100"F,  200”F,  300”F  .  .  . 

Cp(u)  =  20.0,  20.1,  20.2  .  .  . 

4.3.8  Temperature  at  x  =  v 

The  temperature  of  the  compartment  atmosphere  is  a  known  constant  equal  to  u(v). 
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4.3.9  Cross-sectional  Area 

Every  member  contributing  to  the  internal  flow  of  heat  that  does  not  have  the  same 
cross-sectional  area  as  the  representative  element  chosen  is  given  as  Sj. 

4.3.10  Length  of  Heat  Transfer  Path 

The  thickness  or  length  of  every  member  in  the  direction  of  heat  flow  is  given  as  /j. 

4.3.11  Coolant  Enthalpy  Change 

The  heat  sink  capability  of  the  expendable  coolant  is  a  function  of  the  ambient  pressure 
at  which  the  conversion  from  liquid  to  a  gas  takes  place.  This  pressure  can  generally  be 
assumed  to  be  equal  to  that  existing  at  flight  altitudes;  therefore,  it  can  be  related  to 
flight  time.  Example: 

For  0  <  t  <  150  minuses: 

Time  increment.  At  =  150  minutes 

Total  time,  t^  =  150  minutes 

Enthalpy,  hjg  =  1050 

For  150  <  t  <  200: 

Time  increment.  At  =  5  minutes 
Total  time,  tj.  =  50  minutes 
Enthalpy,  h^g  =  1049,  1047  .  .  . 

4.4  Program  Output 

This  section  lists  the  desired  information  resulting  from  an  IBM  7090  Digital 
Computer  solution  of  the  problem  that  is  identified  in  Section  3.4.2. 


4.4.1  Time  Mesh  Size 

The  size  of  the  time  mesh  required  to  satisfy  stability  criteria  must  be  determined. 


4.4.2  Thermal  Property  Coefficients 


The  coefficients  of  the  expressions  for  thermal  conductivity  (Aqj,  and 

specific  heat  (Bqj,  Bj^^)  for  each  material  are  input  requirements. 
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4.4.3  Temperature  Printout 

The  temperature  at  spatial  mesh  points  s,f,  o,  r.andc  atm  times  and  the  temperature 
at  all  spatial  mesh  points  at  the  times  indicated  will  be  of  interest;  for  example  (t^  t^^ 

*30  *100  *150  *m^- 

4.4.4  Heat  Removal  Rates 

The  instantaneous  rates  of  heat  removed  by  the  re-entry  and  compartment  environmental 
cooling  systems  and  at  m  times  are  needed  to  determine  the  and  Q^. 

4.4.5  Total  Heat  Removed 

The  total  amount  of  heat  removed  by  the  re-entry  and  compartment  environmental 
cooling  systems,  and  Q^,  are  required  to  determine  the  coolant  weights. 

4.4.6  Total  Amounts  of  Expended  Coolant 

The  total  weights  of  coolant  expended  by  the  re-entry  and  compartment  environmental 
cooling  systems,  and  W^,  are  used  to  calculate  total  system  weights. 


4.5  Problem  Solution  Procedure 

This  section  describes  a  routine  for  determining  the  solution  of  the  problem  statement 
outlined  in  Section  3.4.2.  If  the  operations  are  performed  in  the  sequence  listed,  the 
temperature  distributions,  rates  of  heat  transfer,  total  heat  removed,  and  the  weight  of 
the  thermal  protection  system  can  be  calculated.  Basically,  the  procedure  is  divided 
into  two  sections;  first,  a  sub-routine,  which  generates  preliminary  information  such  as 
expressions  for  k(u)  and  c  (u),  the  time  mesh  size,  p,  and  the  transformation  of  certain 

input  data  that  are  expressed  as  a  function  of  flight  time  so  that  it  coincides  with  the 
calculated  time  mesh  size;  and  second,  the  order  of  solution  of  the  finite  difference 
equations  derived  in  Sections  3.3  and  3.4.  The  step-by-step  procedures  are  given  in 
the  subparagraphs  that  follow. 


by  the  method  of  least  squares  as  presented  in  Appendix  II. 
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4.5.2  Time  Mesh  Size 

4.5.2. 1  From  the  input  of  surface  temperatures  as  a  function  of  flight  time,  pick  out  the 
maximum  Ug. 

4. 5. 2. 2  Solve  for  the  maximum  k  hy  inserting  Ug  maximum  in  place  of  u  in  the  expression 

2 


“l  =  Aot+  *11“  +  *21“' 


(125) 


4. 5. 2. 3  Solve  for  the  maximum  Cp  by  inserting  Ug  maximum  in  place  of  u  in  the  expression 

(126) 


•^P,  =  ®0I  “ 


4. 5. 2. 4  From  the  input  of  material  densities,  select  p,  and  insert  together  with  h,  kj,  and 


Cpj^  into 


,  Cp, 


2  k. 


(127) 


4. 5. 2. 5  Repeat  steps  in  paragraphs  4. 5. 2. 2  through  4. 5. 2. 4  for  j  =  2  to  determine  P2. 

4. 5. 2. 6  Select  the  smaller  value,  either  p^  or  P2,  and  convert  it  to  the  nearest  preceding 
integer,  which  is  an  exact  divisor  of  the  total  flight  time. 

4. 5. 2. 7  Transform  the  initial  set  of  surface  temperatures,  u(s,t)  into  a  new  set  with  a 
time  increment  equal  to  p  by  the  process  of  interpolation. 

4. 5. 2. 8  Repeat  step  in  Section  4. 5. 2. 7  for  u(r, t)  and  h^g(r,t). 

4.5.3  Solution  of  Finite  Difference  Equations 

4. 5. 3.1  Determine  the  temperatures  for  all  interior  mesh  points  ^2^.  .  .  ., 

^12'  ^22'  ■  ■  ■’  *^n2’  ^13’  ^23'  '  '  ’’  ^n3  m  =  1  from  the  following  equation: 


jm 


P  r 

♦  1  =Uiim+ - 5 - M*ii-''2A,i  Ujim  ) 


*^i+i,i,m  ^^i  +  i,j,in  ^  i-i,  j.m  ^  i-i,i,  m  >  +  4  ( j  +  A,j  U|jm+ 

*2j  ^ijm*  *^i+i,j,m  (128) 


*  Temperatures  at  specific  locations  are  denoted  by  U  in  this  section. 
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4. 5. 3. 2  Determine  the  temperature  at  the  exterior  mesh  points  x  =  f  and  x  =  o  at  time 
m  =  1  from  the  following  equations: 


.C,+  C3+ 


where 


C.  =  (4Ao4+2A,^U,  +2A,^  U,„,+  Aj^  U* 

2''24  U,  U*„) 

C,  =  A,+  ^  U, -H  U,.,  + 


'  4Jt. 


(4Ao,+  2A„  ..m+2A,,  U„^4A„  + 


2A*,  U,  ,„U„  „  -I-  A„ 


Uo.m  )  . 


^4=  Po  (Uf%+  '8  38  .8  +  12  67.9  X  lo’  uf  ,^+  38  8  5 .6  X  1 0*  - 

Uotm-  1838  8  1267.9  X  lo’  U*  „,  -  3885.6  X  hO®  ,  (^34) 


C9  = 


2  ~2~  '^©♦1,  m  +  i  +  ^  ^o,m^  2 


A  2  2 

4  ^  o  ♦  I  ,  m  ♦  I  (1' 

4. 5. 3. 3  Determine  the  rate  of  heat  transfer  at  time  m  =  1  to  and  from  and  the  rate  of 
heat  storage  at  the  exterior  mesh  point  x  =  r.  Now  the  rate  of  heat  transfer  to  r 

<’i9  -4A„U,^„,,  +2A,  U  -2A,U,%,,  +A^U®.,,  + 

^2  ^r-  I,  m*i  t  -  I ,  m  ♦  I  ^r,  m*  I  ~  *2  ^r,  nn»i  )  •  (1 

**0  9  ■  ^0  ^r,m*i  ~^*o^r*i,m»i  ■*‘2A,  2A,  ,  +  Aj  ^  + 

*^r,nn*i'^r  +  i  ,m*  I  ~  ^r,  tn*i^r  ♦! ,  m*i  “A2U^^,^^^,  ^ 
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(138) 


where  U^,  m  +  1  is  given. 

4.5.3.4  Evaluate  conditions  at  x  =  r  and  time  m  +  1.  Three  possibilities  exist.  If 

“i.  •  '»o«+  *«..  •  (139) 

then  no  heat  is  removed  by  the  cooling  system,  or  =  0,  and  U^,  m  + 1  =  U^,  m  + 1 
(given  by  Input). 


4.  ^4  +4 

^05  •» 


b.  When 

then  the  heat  removed  by  the  cooling  system  is 

4 


'ds 


IS  OS  SB 


(140) 

(141) 


and  U^,  m  +  1  =•  U^.  m  +  1  (given  by  input). 


c.  Finally,  if 


q.  <  q  Q 

^IS  ^OS  ^  SB 


(142) 


then  no  heat  is  removed  by  the  cooling  system  (qy,  =  0)  and  U  ,m  +  l^U  ,m  +  l  (given 

i)  r  r 

by  input).  When  condition  (c)  applies  a  new  U^,  m  +  1  must  be  determined  from  the 
following  relationship,  converting  Equation  (142)  to  an  equality 


At 

r.^ 

{‘^r. 

i,m*t  i,m»  i  ] 

+  f- 

2h 

1-  2h 

*  V 

2A« 

A  1 

At 

u‘. 

-  U  * 

)  U 

h 

h 

•4h 

,  m*i  4h  r»i 

,  m*i  /  “r. 

m  ♦  1 

2h  r-i,m*i 


A «  9  A«\  A  9 

4-  ^r-  i,m*i  +  I  ,m  ♦  I  + 


'r  4^  I ,  m  ♦  I 


m  )  =0 

f  ,  01  / 


(143) 


4. 5. 3. 5  Determine  the  temperature  at  the  external  mesh  point  x  =  c  for  m  +  1  from  the 
following  equation: 

^c,m  +  l  ■  +”c  )  ’  (144) 


where 
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4.6  Program  Flow  Chart 

The  program  flow  chart  presents  in  a  concise  manner  a  step-by-step  solution 
procedure  for  a  complex  problem  involving  a  number  of  operations.  Two  flow  charts 
outlining  the  manner  of  solving  the  relationships  for  the  case  in  Section  3.4.2  are  shown 
in  Figures  16  and  17.  When  the  path  connecting  the  events  is  traced,  the  sequence  of 
introducing  or  storing  terms  is  indicated  as  well  as  the  origin  of  each,  input.  A  clearer 
picture  of  the  method  of  solution  may  be  obtained  if  the  flow  chart  of  the  program  is 
used  in  conjunction  with  the  mesh  network. 
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At  this  poiif  some  comments  regarding  the  solution  of  certain  of  the  equations  in 
Section  4.5  should  be  in  order.  In  general,  the  majority  of  the  equations  involve  the 
straightforward  application  of  the  principles  of  algebra  that  result  in  a  rather  obvious 
solution.  On  the  other  hand,  the  solution  of  a  polynomial,  such  as  Equation  (143),  involves 
some  knowledge  of  the  nature  of  the  one  valuable  root.  Of  course,  the  root  of  interest 
must  be  real;  however,  this  condition  by  itself  is  certainly  not  a  sufficient  one  since  it 
is  entirely  feasible  that  the  solution  would  yield  more  than  one  real  root.  Therefore,  a 
second  condition  should  be  included  to  further  restrict  the  number  of  possible  results.  A 
possibility  here  is  to,  in  some  fashion,  bound  the  solution.  For  instance,  if  some  insight 
into  a  limited  range  for  the  one  valuable  root  can  be  applied  by  confining  the  neighborhood 
of  the  solution  within  finite  limits,  the  probability  that  more  than  one  real  root  will  lie 
within  these  bounds  is  low.  In  fact,  the  extremes  usually  can  be  defined  so  that  one  valuable 
root  will  almost  always  be  the  only  one  that  satisfies  both  of  the  previously  mentioned 
restrictions. 


PROGRAM 

INPUT 


Figure  16.  Program  Flow  Chart  for  Subroutine 
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Figure  17.  Program  Flow  Chart  for  the  Insulated  and  Cooled  Compartment 
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5.0  RESULTS 


5.1  System  Thermal  Analysis 

5.1.1  General- 

Numerical  solutions  for  each  of  the  three  problems  presented  in  Section  3  would 
represent  a  considerable  effort  in  view  of  the  number  of  parameters  involved.  Such  a 
treatment  of  each  case  is  not  within  the  scope  of  this  section;  instead,  the  application  of 
the  finite  difference  equations  and  the  method  of  optimizing  the  design  of  the  thermal 
protection  system  will  be  illustrated  for  only  one  case. 

The  example  selected  was  that  designated  in  Section  3.4.1  “The  Insulated  and  Cooled 
Structure.”  This  problem  is  considered  a  compromise  between  the  lengthy  and,  possibly, 
less  meaningful  solutions  that  are  required  for  the  case  in  Section  3.4.2,  and  the  more 
simplified,  straightforward  approach  applying  to  the  case  in  Section  3.4.3.  To  clarify  the 
results,  we  modified  the  case  in  Section  3.4.1  by  terminating  the  solution  at  x  =  c  (see 
Figure  11).  In  addition,  a  thermal-conductance  term  that  depended  on  temperature  was 
introduced  to  include  the  effects  of  an  insulation  package  on  the  heat  transmitted  to  the 
cooling  system.  As  a  result.  Equation  (28)  was  modified  to  incorporate  an  additional 
term. 

5.1.2  Program  Inputs 

This  section  lists  most  of  the  program  inputs  in  tabular  form  as  suggested  in  Section 
4.3.  In  the  case  of  surface  temperatures  of  the  vehicle  during  re-entry,  these  tempera¬ 
tures  can  be  presented  more  conveniently  in  a  graphical  manner  because  of  the  number 
involved.  Optimum  insulation  thicknesses  will  be  determined  for  four  histories  of  surface 
temperatures  U(s,t),  shown  in  Figure  18,  which  are  typical  of  those  experienced  at 
various  locations  on  the  outer  heat  shield.  The  number  of  curves  and  the  maximum 
temperature  of  each  were  selected  to  enable  a  parametric  survey  and  to  facilitate  the 
presentation  of  data.  A  listing  of  the  IBM  program  input  data  is  available  in  Table  1. 

This  table  includes  the  thermal  property  data  determined  from  the  physical  makeup 
of  the  thermal  protection  system.  Coefficients  for  all  polynomials  were  then  determined 
by  the  method  in  Appendix  11.  The  origin  of  most  of  this  information  is  included  in  the 
references. 

Three  thicknesses  of  insulation  were  selected  as  a  minimum  for  demonstrating  the 
influence  of  this  parameter  on  the  weight  of  the  system.  The  smallest  thickness, 
arbitrarily  chosen  as  0.2  inch,  was  established  in  view  of  the  magnitude  of  the  spatial 
mesh  size.  The  time  mesh  size  was  calculated  by  employing  the  stability  criteria 
discussed  in  Section  4.5.2.  Once  this  was  established,  the  surface-temperature  histories, 
through  could  be  tabulated  by  listing  the  temperatures  that  corresponded  to 

each  time  step. 

This  completed  the  input  requirements  for  solving  the  problem. 
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Figure  18.  History  of  Surface-Equilibrium  Temperature  of  Vehicle 
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TABLE  1 


PROGRAM  INPUTS 


Inputs 

Insulation  Thickness 
(inches) 

0.6 

0.8 

1.0 

Number  of  Internal  Mesh  Points  (i) 

6 

7 

8 

Coefficients  of  the  Polynomial  k(u) 

^0 

0.14 

0.14 

0.14 

(Insulation): 

^1 

-  0.177  X  lO'^ 

-  0.177  X  10‘^ 

-  0.177  X  10’^ 

^2 

0.103  X  10'^ 

0.103  X  10'^ 

0.103  X  10’^ 

Coefficients  ofthe  PolynomialCp(u) 

0.23 

0.23 

0.23 

«1 

0 

0 

0 

Insulation  Density  (p) 

12 

12 

12 

Air  Space  Conductance  C  (u) 

^Oa 

1 

0.45 

0.45 

0.45 

A, 

la 

1 

0.4133  X  10'^ 

0.4133  X  10'^ 

0.4133  X  lO'^ 

Heat  Shield  Conductance  C.  (u) 

'^0 

15.642 

15.642 

15.642 

-  19.5  X  lO"^ 

-  19.5  X  lO"^ 

-19.5  X  10'^ 

12.61  X  10‘^ 

12.61  X  10“^ 

12.61  X  lO'^ 

Insulation  Package  Conductance 

5.33  X  lO'^ 

4.0  X  lO"^ 

3.2  X  10‘^ 

Cf  (u) 

^1 

8.667  X  lO"^ 

6.5  X  10'^ 

5.2  X  lO"'’ 

0 

0 

0 

a  F  e 

2.884  X  10"^° 

2.884  X  10'^° 

2.884  X  10”^° 

Spatial  Mesh  Size  (h) 

0.2 

0.2 

0.2 

Time  Mesh  Size  (p) 

0.5 

0.5 

0.5 

Number  of  Time  Steps  (m) 

161 

161 

161 

Initial  Temperatures 

70 

70 

70 

(except  at  X  =  c) 

Temperature  at  x  =  c  (U^q,  .  . 

•) 

200 

200 

200 
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5.1.3  Program  Output 


Numerical  solutions  to  the  problem  were  obtained  by  following  a  procedure  similar 
to  the  one  outlined  in  Section  4.5,  which  led  to  many  of  the  outputs  suggested  in  Section 
4.4. 

The  first  results  are  graphical  representations  of  the  temperature  distribution  through 
two  thicknesses  of  insulation.  0.6  and  1.0  inch,  at  the  re-entry  times  specified.  The 
surface  temperature  history,  which  applies  in  each  case,  has  been  noted  and  the  curves 
are  presented  in  Figures  19  through  26. 
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Figure  19.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 
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Figure  20.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 
of  0.6  Inch  (Surface  Temperature 
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Figure  21.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 
of  0.6  Inch  (Surface  Temperature 
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Figure  22.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 
of  0.6  Inch  (Surface  Temperature 
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Figure  24.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 

of  1.0  Inch  (Surface  Temperature  U  J 
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Figure  25.  Temperature  Gradient  Through  Thermal  Protection  System  for  Insulation  Thickness 
of  1. 0  Inch  (Surface  Temperature  U  ^ 
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Several  interesting  observations  can  be  made  from  these  results.  For  15  minutes  after 
the  start  of  re-entry,  we  note  that  the  temperature  gradient  is  zero  at  some  point  within 
the  insulation.  On  either  side  of  this  point,  the  gradient  indicates  that  heat  is  conducted 
from  both  the  outside  and  the  inside  faces  of  the  insulation  to  internal  elements.  Once  the 
temperature  of  the  external  surface  begins  to  increase,  heat  is  conducted  through  the 
heat  shield  and  the  edges  of  the  insulation  package  to  the  inside  face  of  the  insulation.  The 
response  of  this  path  is  very  rapid  due  to  the  negligible  heat  capacity  of  each  of  these 
components.  On  the  other  hand,  the  heat  capacity  of  the  insulation  cannot  be  neglected; 
therefore,  the  transfer  gradient  through  the  insulation  lags  the  transfer  gradient  through 
the  heat  shield  supports.  Simultaneous  conduction  along  both  paths  gives  the  appearance 
of  a  heat  generation  term  at  the  inside  face  of  the  insulation. 

As  the  time  from  the  start  of  re-entry  increases,  the  point  of  zero  temperature 
gradient  moves  to  the  right  until  it  passes  through  the  insulation.  At  this  time,  heat  is 
being  transferred  in  one  direction  only. 

The  assumption  that  tiie  temperature  of  the  primary  structure  is  initially  130°F  above 
the  temperature  of  the  Insulation  influences  two  directional  heat  flow.  This  assumption 
does  not  affect  the  temperature  distribution  beyond  a  15-minute  period  for  all  but  one 
case,  and,  in  general,  any  contribution  of  this  assumption  to  the  results  is  rather  insignifi¬ 
cant. 

One  will  also  note  that  the  heat  shield  plays  a  very  small  role  in  attenuating  heat 
transfer  to  the  primary  structure.  As  mentioned  previously,  its  major  function  is  to 
transfer  aerodynamically  generated  heat  back  across  the  boundary  layer.  The  equivalent 
conductance  of  this  honeycomb  panel  is  relatively  high;  hence,  the  temperature  gradient 
across  its  thickness  is  small. 

In  comparison,  the  air  space  between  the  insulation  and  the  primary  structure  offers 
considerable  resistance  to  the  flow  of  heat  as  indicated  by  the  temperature  drop  between 
these  two  surfaces. 

The  rates  of  heat  transfer  to  the  cooling  system  are  shown  in  F  igures  27  through  30  for 
each  thickness  of  insulation  as  a  function  of  re-entry  time.  For  each,  the  coolant  require¬ 
ments  are  computed  from  the  value  of  the  integral  under  each  curve.  The  amounts  of 
heat  absorbed  for  each  computer  run  are  given  in  Table  2. 

TABLE  2 

HEAT  ABSORBED  FOR  COMPUTER  RUNS 

Total  Heat  Absorbed  (BTU  per  Sq  Ft)  ) 


Temperature  History 


Insulation  Thickness 
(inches) 

^S1 

Us2 

^s3 

^S4 

0.6 

454 

355 

149 

44 

0.8 

325 

248 

101 

27 

1.0 

232 

174 

68 

15 

59 


RE-  ENTRY  TIME  (MINUTES) 


Figure  28.  Etlectof  Insulation  Thickness  on  Rate  of  Heat  Transfer  to  Coolant  (Surface  Temperature  U  ^ 
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The  maximum  rate  of  heat  transfer  is  shown  in  Figure  31  and  the  coolant  requirements 
are  shown  in  Figure  32  as  each  is  affected  by  the  insulation  thickness  and  the  maximum 
equilibrium  temperature.  The  use  of  this  latter  term  as  a  parameter  was  justified,  since 
all  temperature  histories  are  of  the  same  duration  and  the  general  shape  of  the  curves 
is  similar. 

Figure  32  indicates  the  heat  transmitted  to  the  coolant  is  nearly  zero  if  the  maximum 
surface  equilibrium  temperature  is  not  high  for  some  insulation  thicknesses.  For  this 
case,  the  elimination  of  the  cooling  system  might  be  considered. 


Figure  31.  Effect  of  Temperature  and  Insulation  Thickness  on  Rate  of  Heat  Transfer 
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Figure  32.  Effect  of  Temperature  and  Insulation  Thickness  on  Coolant  Weight 


5.2  Insulation  and  Coolant  Tradeoff 

If  the  coolant  requirements,  as  shown  in  Figure  32,  are  combined  with  the  weight  of 
the  insulation,  a  relationship  between  the  total  weight  of  the  thermal  protection  system 
and  the  thickness  of  insulation  can  be  demonstrated  for  a  number  of  outer  surface 
temperatures.  These  results  are  shown  in  Figure  33  in  which  the  maximum  surface 
temperature  is  used  as  a  notation  of  the  particular  history  (refer  to  Figure  18).  An 
examination  of  these  results  reveals  that  the  selection  of  a  thickness  of  insulation  can 
be  made,  which,  when  combined  with  the  amount  of  coolant  required,  will  result  in  a 
thermal  protection  system  of  minimum  weight. 

This  type  of  presentation  may  be  conveniently  used  as  design  information,  since  the 
surface  of  the  re-entry  vehicle  can  be  sectioned  into  temperature  zones.  The  boundaries 
of  these  zones  can  be  judiciously  assigned,  based  on  temperature  distributions  and  the 
surface  temperature  history  for  the  re-entry  trajectory,  which  results  in  maximum 
heating.  In  the  case  of  glide  vehicles,  temperature  gradients  on  the  surface  are  relatively 
insignificant;  hence,  the  total  number  of  zones  would  be  few  and  can  be  based  on  a  modest 
deviation  from  an  average  temperature.  As  a  result,  a  procedure  is  established  that  finds 
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application  as  a  practical  design  tool.  The  more  problematical  aspects  are  those  that  are 
associated  with  the  derivation  of  analytical  relationships  and  the  programming  of  these  in 
;ach  case.  These  areas  were  emphasized  in  this  study. 


Figure  33.  Insulation  -  Coolant  Weight  Tradeoff 
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APPENDIX  I 

DERIVATION  OF  FINITE  DIFFERENCE  RELATIONSHIPS 


In  solving  a  differential  equation  of  the  form  in  Equation  (13),  one  can  express 
conveniently  the  temperature  distribution,  U^.  as  an  nth  degree  polynomial  containing 
(n  +  1)  arbitrary  constants  as 

(I)  U)  (3)  (n) 

IJk  •=  “o  +  “i  *  +  02  «  -H  Q,  «  + .  +  “n  *  (152) 

If  X  is  replaced  by  the  variable 

*  —  xo 

“  "  (153) 

where  h  is  the  interval  between  two  successive  values  of  the  argument  x,  the  differences 

of  U,,  may  now  be  used  to  determine  the  values  of  the  coefficients  a^,,  a, ,  a„,  .  .  .  ,  a 
V  ■'  U  1  2  n 

thus 

.  ('1  (n-i) 

AU^  sai+aojv  +3a3w  +  --+nqn 

2  O)  (n-2) 

A  Uy  =  2- !•  Oj  +  3-2- a  3  V  +■•..  +n(n-i)an 

A’u^  =3-21-a  + -  +  o|n-i)(n-2) 


a" Uy  =  a  „  (  n  /  )  • 


Since  these  expressions  must  hold  for  all  values  of  x,  the  coefficients  can  now  be  found 
by  setting  v  =  o.  Therefore 


=  U, 


=  AU, 


2  / 


A  Up 
n  / 


Substituting  the  values  for  the  coefficients  into  Equation  (152),  Gregory-Newton’s 
Interpolation  formula  is  obtained  thus 


Uv  =  'Jo  + 


0) 


AUo  + 


2  / 


A*  U„ 


n  /  o 


(155) 


Now  the  value  of  a  function  of  x  may  be  found  by  replacing  the  function  with  this  inter¬ 
polating  polynomial,  then  differentiating 


d  f ( X )  ^  d  U  X  dUy  dv 

dx  dx  dv  dx  X  =  X  (1^^) 
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The  derivative  of  Equation  (153)  gives 

dw  _  _l_ 
d  X  h 

Substituting  Equation  (157)  into  Equation  (156)  gives  the  derivative 


(157) 


df  (  x)  ^  I  dUy 

dx  “  h  dv  (158) 


An  approximation  of  the  second  derivative  of  f  (x)  can  be  found  by  again  differentiating 
Equation  (156) 


d®  f  (  X  ) 

_  d"*  U,  dpi  dU„  - 

1  d  v 

dx  * 

“2  L 

dx  dv  h  dv  ' 

'  dx 

d  V  1 

Recalling  that  — ^ - =  -jj — 

,  Equation  (159)  can  be  rewritten 

(159) 


d ‘  f  (5) 


d  x 


d  V 


(160) 


Applying  the  expression  relating  the  derivative  of  the  function  to  the  derivative  of  the 
interpolation  formula  (155)  gives  the  first  derivative 


d  u  ^  I  r  d  U  V  1  _ 

dx  h  ^  dv  v=o  ^ 


+  —  A^Uc 
3  ‘ 


- A  u. 


Clearly,  the  first  approximation  is 


or  in  general 


X=  I 


’•(161) 


(162) 


This  approximation  is  most  usefully  applied  to  the  boundary  condition  where  the 
derivative  does  not  exist  for  the  i-1  element. 


The  finite  difference  approximation  for  the  time  derivative  is  found  in  a  similar 
manner  resulting  in 


(tt)  « -t  I  ■ 


t=t. 


Here  p  is  the  difference  between  two  successive  values  of  the  argument  t. 


(163) 
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Another  interpolation  formula  can  be  developed  by  averaging  the  forward  and  backward 
interpolation  formulas  of  Gauss.  The  forward  interpolation  formula  expresses  the 
temperature 


u  =  Uo+  vAU,/^  +  ^  (  V  +  I)‘®’  a’  + 


(  zn  -  I)  / 


(z  n  - 1  )  2  n-i 

(  V  +  n  -I  )  A  U 


'/z 


I 


( z  n  1  2  n 

(  V  4-  n-  I )  A  U_  + 


(  z  n)/ 

In  the  backward  interpolation  formula 


(164) 


u  =  U-  +  vAU.i,  +  -i-  (  w  +  1)'“  A‘  U-  +  — ,  (v  +  D'^'A'’  U.y,  + 
o  ^2  2/  O  /  '2 


A* 


_l_ 

3/ 


,(3)*5 


I 


(2n)/ 


«  A*”  II 

(w  -»-n  )  A  Uft 


t 


(2n+  I)  / 


(  V  +  n  ) 


(z  n  +  I  ) 


an  ♦  1 

A  4- 


(165) 


Averaging  Equations  (164)  and  (165)  gives 


AU.,  +AU„ 


V*  A* 

2  *  “  --I 


I  z  A*  U  _  2  ^  1 

U  ,  +  4-  V  (  V  -I)  - =-^ - =J-  + 

6  2 


(166) 


which  is  the  central  difference  formula  of  Stirling. 

If  the  relationship  (158)  is  applied  to  Equation  (166),  an  approximation  of  the  first 
derivative  may  be  determined 


/  du  \  „  1  r  <lUv  1  (  /  +AUo 

1  a*u_2  +  A*u., 

.V 

'dx'  “hL(|v  J“h  V  2 

*-*o  v=o 

6  2 

■J  (167) 

Again,  using  the  first  approximation  gives 

,du  .  .. 

U,  -U., 

(168) 

(dx  )  .  "  2h 

2h 

Now  the  derivative  squared  term  is  obtained  simply  by  squaring  the  right-hand  member 
of  Equation  (168) 


_J _ 

4h^ 


(U,  -2U, 


(169) 
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From  Equation  (166),  the  first  approximation  of  the  second  derivative  is 

9  2 

/  d  u  \  ^  r  ^  1  I  2 

)  ~  -r  - ^  1  =  ,  ) 

'  d  *  2  Xr  ^  ^  Hw®  J..  _  ^ 


{  AU.,)  =  ^  [(  U,  -U„)-(U^  -U.,)  ] 


^  dx*  i 


<  U,  -  2Uo  +  U.,  ) 


(170) 


This  completes  the  derivation  of  all  finite  difference  equations  that  are  required  to  solve 
the  general  heat-transfer  equation  and  the  boundary  conditions  in  Section  3. 


If  the  heat  diffusion  Equation  (13)  could  be  solved  conveniently  by  a  Fourier  series 
expansion,  we  could  show  that  stability  of  the  solution  depends  on  the  ratio 


M 


g  P 
h* 


(171) 


Since  this  would  be  somewhat  lengthy  to  show,  it  will  just  be  noted  that  an  examination 
of  this  parameter  is  required  prior  to  attempting  a  solution  of  the  problem.  If,  for 
instance,  M  <  1/2,  each  term  in  the  solution  of  the  difference  equation  will  decrease 
exponentially  with  time.  This  means  that  round  off  error,  although  carried  to  subsequent 
rows,  would  eventually  vanish.  On  the  other  hand,  if  M  >  1/2,  some  terms  would  increase 
exponentially  with  time,  and  round  off  errors  would  tend  to  be  amplified  as  the  solution 
proceeds.  Therefore,  once  an  h  has  been  selected,  it  becomes  important  to  choose  a 
value  for  p  such  that  M  <  1/2  to  insure  stability  of  the  solution. 
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APPENDIX  II 

METHOD  OF  LEAST  SQUARES 


This  section  reviews  a  method  of  determining  the  constants  that  appear  in  the  equation 
that  is  selected  to  represent  thermal  conductivity  and  specific  heat  data.  It  is  the  most 
useful  method  and  the  one  most  frequently  applied  since  it  has  the  advantage  of  producing 
a  unique  set  of  values  for  the  constants  of  the  polynomials.  Moreover,  these  constants 
give  the  most  probable  equation  in  the  sense  that  the  computed  values  of  k(u)  or  Cp(u) 
are  the  most  probable  values  of  tbe  observations,  since  the  residuals  are  assumed  to 
follow  the  Gaussian  law  of  error.  In  other  words,  the  principle  of  the  method  of  least 
squares  asserts  that  the  most  representative  curve  is  that  for  which  the  sum  of  the 
squares  of  the  residuals  is  a  minimum. 

Suppose  the  given  set  of  observed  values  (k.,  u.),  (i  =  1,  2,  3,  .  .  .  n),  can  be  represented 
by  the  equation 


k  =  f  (  u  ) 

th 

containing  r  undetermined  constants,  Aj,  ....  A^.  Then  the  n  observation 
equations,  for  example 


(172) 


are  to  be  solved  for  the  r  unknowns.  If  r  =  n  there  are  just  enough  conditions  to  determine 
the  constants;  if  n  <  r,  there  are  not  enough  conditions  and  the  problem  is  indeterminate; 
but,  in  general,  n>r,  and  there  are  more  conditions  than  there  are  unknowns.  In  the 
general  case,  the  values  of  an  m  which  satisfy  any  r  of  these  equations  will  not  satisfy 
the  remaining  n  -  r  equations,  and  the  problem  is  to  determine  the  set  of  values  of 
that  will  give  the  most  probable  values  of  k.  Let 


=  k:  -  k 


(173) 


be  the  residuals  or  deviations  of  the  computed  values  from  the  observed  values,  where 
kj  is  the  value  of  k  obtained  by  substituting  u  =  uj  in  k  =  f(u).  On  the  basis  of  the  Gaussian 
law  of  error,  the  probability  of  obtaining  the  observed  values  ki  is 


P 


-  h' 


(174) 


P  is  a  maximum  where  X  v,  is  a  minimum.  Since  s  =  V|  is  a  function 

i  =  I  i  =  I 

of  the  r  unknowns,  A^,  A2.  .  .  .  ,  A^.,  it  follows  that  the  necessary  conditions  for  a 

minimum  are 


d  a 
<3  A, 


=  O  . 


(175) 
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Moreover,  each  Vj  is  a  function  of  a^^;  therefore 


O  ,2.2.  ,  2  , 

‘  ^  . 


d  s  d  V, 

-  =  2v,  -- — — 

d  A,  '  dAm 


2v,  - ^  + 

*  dAr,, 


2  I  v. 


,  (  m  =  1,2,  •  •  •  ,  r  ) 


‘I  1=1 


These  formulas.  Equation  (175).  are  called  the  normal  equation. 

If  the  r  functions  in  Equation  (175)  are  linear  in  the  r  unknown  Aj,  A2. . .  ..A^.then 

these  equations  can  be  solved  immediately.  This  will  certainly  be  the  case  if  f(u)  is  a 
polynomial.  Let 


so  that 


f  (o)  =  Z  U 

i  =  ' 


Vi  =  I  AjUj  -k 

j  =  ' 


Then  —  =  u,  ,  and  the  normal  equations  assume  the  form,  with  the  aid  of 


Equation  (176), 


I  (  i  ''iu,'".,, )  .0, 


(m  =1,  2  ,•••  r  )  • 


1=1  1=1 


One  should  note  that  the  equation  that  is  obtained  by  setting  m  =  1  is 

n 

I  Vj  =  O 
i  =  I 

Rearranging  the  terms  in  Equation  (180)  and  collecting  the  coefficients  of  A^  gives 


'  ,  ”  j  +  m  -  2 


_2:(lu,  )Aj=.ZUj' 


(  m  =  I,  2  ,  ••  •,  r  ) 
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Now  the  r  linear  equations  can  be  solved  for  the  values  of  the  r  unknowns  Aj^,  A^, 


.  .  .  .  A^.  For  example: 


U 


K 


1 

2 

3 

4 


1.7 

1.8 
2.3 
3.2 


2  2 

Now  express  f(u)  =  +  A2U  +  A^u  .  then  v.  =  A^  +  A2n.  +  A^u.  -  k.,  and 
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The  normal  equations  are 
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or  in  other  terms 
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2  5:  Vj  =  O  . 
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I  V  — ^  =  0,  (  m=  I.  2,3  ) 
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I  (A,+  A2Uj+  A,U.*  -  k,  )  ■  I  =  O  . 
i  =1 


Z  (  A,  +Aj  Uj  +  A,  Uj*  -  kj  1  Uj  =  O  , 
i  =  1 


Z  (A,+A,  Ui+  A,U?  -k.)-  Uf  =0 
i  =  I 

If  the  coefficients  of  A^  are  collected  and  the  normal  equations  put  in  the  form  of 
Equation  (182),  the  following  three  equations  result; 

4A,  4-  (  Z  Uj  )  Ag  +  (  Z  uf  )a3  =  Z  kj  , 

1=1  in  i=  I 

(  iuj)  A,+  (  Zuj)A,+  (Zu|)a,  =  Z  Ui  ki  , 
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Now 


so  that 


(  Zui*)A.+  (  (Iur)A3  =  I  Ui* 

'  I  =1  'in  '  '  1=1  '  1=1 


Zuj  =  I+2+3  +  4  =  10, 

1  =  1 

I  Uj  =  1  +  4+9+16  =  30,  •»c. 

i  =1 


4  A,  +  10  Aj  •+  30A  3  =  9  , 

10  A,  +  3OA3  +IOOA3  =25, 
30A,  +IOOA2  +354A3  *  80.8 
When  these  equations  are  solved  for  the  coefficients. 


A,  =  2  ,  A,  =  -0.  5,  A; 
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=  0.2. 
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